diff --git a/src/DynamicModel.cc b/src/DynamicModel.cc
index 097b8a736d753b7840e83552421a0852b2946da5..6a7b5dd6841cc5c6a5b1b2e3654f01b6696a4d94 100644
--- a/src/DynamicModel.cc
+++ b/src/DynamicModel.cc
@@ -4822,27 +4822,25 @@ DynamicModel::computingPass(bool jacobianExo, int derivsOrder, int paramsDerivsO
first_order_endo_derivatives = collect_first_order_derivatives_endogenous();
is_equation_linear = equationLinear(first_order_endo_derivatives);
- evaluateAndReduceJacobian(eval_context, contemporaneous_jacobian, static_jacobian, dynamic_jacobian, cutoff, false);
+ tie(contemporaneous_jacobian, static_jacobian) = evaluateAndReduceJacobian(eval_context, cutoff, false);
if (!computeNaturalNormalization())
- computeNonSingularNormalization(contemporaneous_jacobian, cutoff, static_jacobian, dynamic_jacobian);
+ computeNonSingularNormalization(contemporaneous_jacobian, cutoff, static_jacobian);
lag_lead_vector_t equation_lag_lead, variable_lag_lead;
- blocks = select_non_linear_equations_and_variables(is_equation_linear, dynamic_jacobian, equation_reordered, variable_reordered,
- inv_equation_reordered, inv_variable_reordered,
- equation_lag_lead, variable_lag_lead,
- n_static, n_forward, n_backward, n_mixed);
+ tie(blocks, equation_lag_lead, variable_lag_lead, n_static, n_forward, n_backward, n_mixed)
+ = select_non_linear_equations_and_variables(is_equation_linear);
- equation_type_and_normalized_equation = equationTypeDetermination(first_order_endo_derivatives, variable_reordered, equation_reordered, 0);
+ equation_type_and_normalized_equation = equationTypeDetermination(first_order_endo_derivatives, 0);
prologue = 0;
epilogue = 0;
- block_type_firstequation_size_mfs = reduceBlocksAndTypeDetermination(dynamic_jacobian, blocks, equation_type_and_normalized_equation, variable_reordered, equation_reordered, n_static, n_forward, n_backward, n_mixed, block_col_type, linear_decomposition);
+ block_type_firstequation_size_mfs = reduceBlocksAndTypeDetermination(blocks, equation_type_and_normalized_equation, n_static, n_forward, n_backward, n_mixed, linear_decomposition);
- computeChainRuleJacobian(blocks_derivatives);
+ computeChainRuleJacobian();
- blocks_linear = BlockLinear(blocks_derivatives, variable_reordered);
+ blocks_linear = BlockLinear();
collect_block_first_order_derivatives();
@@ -4855,29 +4853,29 @@ DynamicModel::computingPass(bool jacobianExo, int derivsOrder, int paramsDerivsO
if (block)
{
- evaluateAndReduceJacobian(eval_context, contemporaneous_jacobian, static_jacobian, dynamic_jacobian, cutoff, false);
+ tie(contemporaneous_jacobian, static_jacobian) = evaluateAndReduceJacobian(eval_context, cutoff, false);
- computeNonSingularNormalization(contemporaneous_jacobian, cutoff, static_jacobian, dynamic_jacobian);
+ computeNonSingularNormalization(contemporaneous_jacobian, cutoff, static_jacobian);
- computePrologueAndEpilogue(static_jacobian, equation_reordered, variable_reordered);
+ computePrologueAndEpilogue(static_jacobian);
first_order_endo_derivatives = collect_first_order_derivatives_endogenous();
- equation_type_and_normalized_equation = equationTypeDetermination(first_order_endo_derivatives, variable_reordered, equation_reordered, mfs);
+ equation_type_and_normalized_equation = equationTypeDetermination(first_order_endo_derivatives, mfs);
cout << "Finding the optimal block decomposition of the model ..." << endl;
lag_lead_vector_t equation_lag_lead, variable_lag_lead;
- computeBlockDecompositionAndFeedbackVariablesForEachBlock(static_jacobian, dynamic_jacobian, equation_reordered, variable_reordered, blocks, equation_type_and_normalized_equation, false, true, mfs, inv_equation_reordered, inv_variable_reordered, equation_lag_lead, variable_lag_lead, n_static, n_forward, n_backward, n_mixed);
+ tie(blocks, equation_lag_lead, variable_lag_lead, n_static, n_forward, n_backward, n_mixed) = computeBlockDecompositionAndFeedbackVariablesForEachBlock(static_jacobian, equation_type_and_normalized_equation, false, true);
- block_type_firstequation_size_mfs = reduceBlocksAndTypeDetermination(dynamic_jacobian, blocks, equation_type_and_normalized_equation, variable_reordered, equation_reordered, n_static, n_forward, n_backward, n_mixed, block_col_type, linear_decomposition);
+ block_type_firstequation_size_mfs = reduceBlocksAndTypeDetermination(blocks, equation_type_and_normalized_equation, n_static, n_forward, n_backward, n_mixed, linear_decomposition);
printBlockDecomposition(blocks);
- computeChainRuleJacobian(blocks_derivatives);
+ computeChainRuleJacobian();
- blocks_linear = BlockLinear(blocks_derivatives, variable_reordered);
+ blocks_linear = BlockLinear();
collect_block_first_order_derivatives();
@@ -4888,7 +4886,8 @@ DynamicModel::computingPass(bool jacobianExo, int derivsOrder, int paramsDerivsO
computeTemporaryTermsOrdered();
int k = 0;
equation_block.resize(equations.size());
- variable_block_lead_lag = vector<tuple<int, int, int>>(equations.size());
+ variable_block_lead_lag.clear();
+ variable_block_lead_lag.resize(equations.size());
for (unsigned int i = 0; i < getNbBlocks(); i++)
{
for (unsigned int j = 0; j < getBlockSize(i); j++)
@@ -5072,11 +5071,11 @@ DynamicModel::get_Derivatives(int block)
}
void
-DynamicModel::computeChainRuleJacobian(blocks_derivatives_t &blocks_endo_derivatives)
+DynamicModel::computeChainRuleJacobian()
{
map<int, expr_t> recursive_variables;
unsigned int nb_blocks = getNbBlocks();
- blocks_endo_derivatives = blocks_derivatives_t(nb_blocks);
+ blocks_derivatives.resize(nb_blocks);
for (unsigned int block = 0; block < nb_blocks; block++)
{
block_derivatives_equation_variable_laglead_nodeid_t tmp_derivatives;
@@ -5084,7 +5083,6 @@ DynamicModel::computeChainRuleJacobian(blocks_derivatives_t &blocks_endo_derivat
int block_size = getBlockSize(block);
int block_nb_mfs = getBlockMfs(block);
int block_nb_recursives = block_size - block_nb_mfs;
- blocks_endo_derivatives.push_back(block_derivatives_equation_variable_laglead_nodeid_t(0));
for (int i = 0; i < block_nb_recursives; i++)
{
if (getBlockEquationType(block, i) == E_EVALUATE_S)
@@ -5110,7 +5108,7 @@ DynamicModel::computeChainRuleJacobian(blocks_derivatives_t &blocks_endo_derivat
}
tmp_derivatives.emplace_back(eq, var, lag, first_chain_rule_derivatives[{ eqr, varr, lag }]);
}
- blocks_endo_derivatives[block] = tmp_derivatives;
+ blocks_derivatives[block] = tmp_derivatives;
}
}
diff --git a/src/DynamicModel.hh b/src/DynamicModel.hh
index 4dacd05852385afd8ff76d47cfe6cd9c75fb49df..6033ab7817c0757040c04d10404fb01c06e5e547 100644
--- a/src/DynamicModel.hh
+++ b/src/DynamicModel.hh
@@ -138,7 +138,7 @@ private:
//! return a map on the block jacobian
map<tuple<int, int, int, int, int>, int> get_Derivatives(int block);
//! Computes chain rule derivatives of the Jacobian w.r. to endogenous variables
- void computeChainRuleJacobian(blocks_derivatives_t &blocks_derivatives);
+ void computeChainRuleJacobian();
string reform(const string &name) const;
diff --git a/src/ModelTree.cc b/src/ModelTree.cc
index ca721cf96ed4b9f649d7581379f41d9c1726e9de..71b8550d104a08c7a70d533fbf76346fc6a25441 100644
--- a/src/ModelTree.cc
+++ b/src/ModelTree.cc
@@ -349,7 +349,7 @@ ModelTree::computeNormalization(const jacob_map_t &contemporaneous_jacobian, boo
}
void
-ModelTree::computeNonSingularNormalization(jacob_map_t &contemporaneous_jacobian, double cutoff, jacob_map_t &static_jacobian, dynamic_jacob_map_t &dynamic_jacobian)
+ModelTree::computeNonSingularNormalization(jacob_map_t &contemporaneous_jacobian, double cutoff, jacob_map_t &static_jacobian)
{
bool check = false;
@@ -466,9 +466,10 @@ ModelTree::computeNormalizedEquations() const
return endo2eqs;
}
-void
-ModelTree::evaluateAndReduceJacobian(const eval_context_t &eval_context, jacob_map_t &contemporaneous_jacobian, jacob_map_t &static_jacobian, dynamic_jacob_map_t &dynamic_jacobian, double cutoff, bool verbose)
+pair<ModelTree::jacob_map_t, ModelTree::jacob_map_t>
+ModelTree::evaluateAndReduceJacobian(const eval_context_t &eval_context, double cutoff, bool verbose)
{
+ jacob_map_t contemporaneous_jacobian, static_jacobian;
int nb_elements_contemparenous_Jacobian = 0;
set<vector<int>> jacobian_elements_to_delete;
for (const auto &[indices, d1] : derivatives[1])
@@ -527,22 +528,20 @@ ModelTree::evaluateAndReduceJacobian(const eval_context_t &eval_context, jacob_m
cout << jacobian_elements_to_delete.size() << " elements among " << derivatives[1].size() << " in the incidence matrices are below the cutoff (" << cutoff << ") and are discarded" << endl
<< "The contemporaneous incidence matrix has " << nb_elements_contemparenous_Jacobian << " elements" << endl;
}
+
+ return { contemporaneous_jacobian, static_jacobian };
}
-vector<pair<int, int>>
-ModelTree::select_non_linear_equations_and_variables(vector<bool> is_equation_linear, const dynamic_jacob_map_t &dynamic_jacobian, vector<int> &equation_reordered, vector<int> &variable_reordered,
- vector<int> &inv_equation_reordered, vector<int> &inv_variable_reordered,
- lag_lead_vector_t &equation_lag_lead, lag_lead_vector_t &variable_lag_lead,
- vector<unsigned int> &n_static, vector<unsigned int> &n_forward, vector<unsigned int> &n_backward, vector<unsigned int> &n_mixed)
+tuple<vector<pair<int, int>>, lag_lead_vector_t, lag_lead_vector_t,
+ vector<unsigned int>, vector<unsigned int>, vector<unsigned int>, vector<unsigned int>>
+ModelTree::select_non_linear_equations_and_variables(const vector<bool> &is_equation_linear)
{
vector<int> eq2endo(equations.size(), 0);
- /*equation_reordered.resize(equations.size());
- variable_reordered.resize(equations.size());*/
unsigned int num = 0;
for (auto it : endo2eq)
if (!is_equation_linear[it])
num++;
- vector<int> endo2block = vector<int>(endo2eq.size(), 1);
+ vector<int> endo2block(endo2eq.size(), 1);
vector<pair<set<int>, pair<set<int>, vector<int>>>> components_set(num);
int i = 0, j = 0;
for (auto it : endo2eq)
@@ -555,11 +554,9 @@ ModelTree::select_non_linear_equations_and_variables(vector<bool> is_equation_li
i++;
j++;
}
- getVariableLeadLagByBlock(dynamic_jacobian, endo2block, endo2block.size(), equation_lag_lead, variable_lag_lead, equation_reordered, variable_reordered);
- n_static = vector<unsigned int>(endo2eq.size(), 0);
- n_forward = vector<unsigned int>(endo2eq.size(), 0);
- n_backward = vector<unsigned int>(endo2eq.size(), 0);
- n_mixed = vector<unsigned int>(endo2eq.size(), 0);
+ auto [equation_lag_lead, variable_lag_lead] = getVariableLeadLagByBlock(endo2block, endo2block.size());
+ vector<unsigned int> n_static(endo2eq.size(), 0), n_forward(endo2eq.size(), 0),
+ n_backward(endo2eq.size(), 0), n_mixed(endo2eq.size(), 0);
for (unsigned int i = 0; i < endo2eq.size(); i++)
{
if (variable_lag_lead[variable_reordered[i]].first != 0 && variable_lag_lead[variable_reordered[i]].second != 0)
@@ -581,7 +578,8 @@ ModelTree::select_non_linear_equations_and_variables(vector<bool> is_equation_li
inv_variable_reordered[variable_reordered[i]] = i;
inv_equation_reordered[equation_reordered[i]] = i;
}
- return blocks;
+ return { blocks, equation_lag_lead, variable_lag_lead,
+ n_static, n_forward, n_backward, n_mixed };
}
bool
@@ -613,7 +611,7 @@ ModelTree::computeNaturalNormalization()
}
void
-ModelTree::computePrologueAndEpilogue(const jacob_map_t &static_jacobian_arg, vector<int> &equation_reordered, vector<int> &variable_reordered)
+ModelTree::computePrologueAndEpilogue(const jacob_map_t &static_jacobian_arg)
{
vector<int> eq2endo(equations.size(), 0);
equation_reordered.resize(equations.size());
@@ -731,7 +729,7 @@ ModelTree::computePrologueAndEpilogue(const jacob_map_t &static_jacobian_arg, ve
}
equation_type_and_normalized_equation_t
-ModelTree::equationTypeDetermination(const map<tuple<int, int, int>, expr_t> &first_order_endo_derivatives, const vector<int> &Index_Var_IM, const vector<int> &Index_Equ_IM, int mfs) const
+ModelTree::equationTypeDetermination(const map<tuple<int, int, int>, expr_t> &first_order_endo_derivatives, int mfs) const
{
expr_t lhs;
BinaryOpNode *eq_node;
@@ -739,20 +737,20 @@ ModelTree::equationTypeDetermination(const map<tuple<int, int, int>, expr_t> &fi
equation_type_and_normalized_equation_t V_Equation_Simulation_Type(equations.size());
for (unsigned int i = 0; i < equations.size(); i++)
{
- int eq = Index_Equ_IM[i];
- int var = Index_Var_IM[i];
+ int eq = equation_reordered[i];
+ int var = variable_reordered[i];
eq_node = equations[eq];
lhs = eq_node->arg1;
Equation_Simulation_Type = E_SOLVE;
- auto derivative = first_order_endo_derivatives.find({ eq, var, 0 });
pair<bool, expr_t> res;
- if (derivative != first_order_endo_derivatives.end())
+ if (auto derivative = first_order_endo_derivatives.find({ eq, var, 0 });
+ derivative != first_order_endo_derivatives.end())
{
set<pair<int, int>> result;
derivative->second->collectEndogenous(result);
auto d_endo_variable = result.find({ var, 0 });
//Determine whether the equation could be evaluated rather than to be solved
- if (lhs->isVariableNodeEqualTo(SymbolType::endogenous, Index_Var_IM[i], 0) && derivative->second->isNumConstNodeEqualTo(1))
+ if (lhs->isVariableNodeEqualTo(SymbolType::endogenous, variable_reordered[i], 0) && derivative->second->isNumConstNodeEqualTo(1))
Equation_Simulation_Type = E_EVALUATE;
else
{
@@ -775,12 +773,11 @@ ModelTree::equationTypeDetermination(const map<tuple<int, int, int>, expr_t> &fi
return V_Equation_Simulation_Type;
}
-void
-ModelTree::getVariableLeadLagByBlock(const dynamic_jacob_map_t &dynamic_jacobian, const vector<int> &components_set, int nb_blck_sim, lag_lead_vector_t &equation_lead_lag, lag_lead_vector_t &variable_lead_lag, const vector<int> &equation_reordered, const vector<int> &variable_reordered) const
+pair<lag_lead_vector_t, lag_lead_vector_t>
+ModelTree::getVariableLeadLagByBlock(const vector<int> &components_set, int nb_blck_sim) const
{
int nb_endo = symbol_table.endo_nbr();
- variable_lead_lag = lag_lead_vector_t(nb_endo, { 0, 0 });
- equation_lead_lag = lag_lead_vector_t(nb_endo, { 0, 0 });
+ lag_lead_vector_t variable_lead_lag(nb_endo, { 0, 0 }), equation_lead_lag(nb_endo, { 0, 0 });
vector<int> variable_blck(nb_endo), equation_blck(nb_endo);
for (int i = 0; i < nb_endo; i++)
{
@@ -815,10 +812,12 @@ ModelTree::getVariableLeadLagByBlock(const dynamic_jacob_map_t &dynamic_jacobian
equation_lead_lag[j_1] = { -lag, equation_lead_lag[j_1].second };
}
}
+ return { equation_lead_lag, variable_lead_lag };
}
-void
-ModelTree::computeBlockDecompositionAndFeedbackVariablesForEachBlock(const jacob_map_t &static_jacobian, const dynamic_jacob_map_t &dynamic_jacobian, vector<int> &equation_reordered, vector<int> &variable_reordered, vector<pair<int, int>> &blocks, const equation_type_and_normalized_equation_t &Equation_Type, bool verbose_, bool select_feedback_variable, int mfs, vector<int> &inv_equation_reordered, vector<int> &inv_variable_reordered, lag_lead_vector_t &equation_lag_lead, lag_lead_vector_t &variable_lag_lead, vector<unsigned int> &n_static, vector<unsigned int> &n_forward, vector<unsigned int> &n_backward, vector<unsigned int> &n_mixed) const
+tuple<vector<pair<int, int>>, lag_lead_vector_t, lag_lead_vector_t,
+ vector<unsigned int>, vector<unsigned int>, vector<unsigned int>, vector<unsigned int>>
+ModelTree::computeBlockDecompositionAndFeedbackVariablesForEachBlock(const jacob_map_t &static_jacobian, const equation_type_and_normalized_equation_t &Equation_Type, bool verbose_, bool select_feedback_variable)
{
int nb_var = variable_reordered.size();
int n = nb_var - prologue - epilogue;
@@ -865,7 +864,7 @@ ModelTree::computeBlockDecompositionAndFeedbackVariablesForEachBlock(const jacob
// Compute strongly connected components
int num = strong_components(G2, endo2block_map);
- blocks = vector<pair<int, int>>(num, { 0, 0 });
+ vector<pair<int, int>> blocks(num, { 0, 0 });
// Create directed acyclic graph associated to the strongly connected components
using DirectedGraph = boost::adjacency_list<boost::vecS, boost::vecS, boost::directedS>;
@@ -905,7 +904,7 @@ ModelTree::computeBlockDecompositionAndFeedbackVariablesForEachBlock(const jacob
get<0>(components_set[endo2block[i]]).insert(i);
}
- getVariableLeadLagByBlock(dynamic_jacobian, endo2block, num, equation_lag_lead, variable_lag_lead, equation_reordered, variable_reordered);
+ auto [equation_lag_lead, variable_lag_lead] = getVariableLeadLagByBlock(endo2block, num);
vector<int> tmp_equation_reordered(equation_reordered), tmp_variable_reordered(variable_reordered);
int order = prologue;
@@ -927,10 +926,8 @@ ModelTree::computeBlockDecompositionAndFeedbackVariablesForEachBlock(const jacob
add_edge(vertex(i, G2), vertex(i, G2), G2);
//Determines the dynamic structure of each equation
- n_static = vector<unsigned int>(prologue+num+epilogue, 0);
- n_forward = vector<unsigned int>(prologue+num+epilogue, 0);
- n_backward = vector<unsigned int>(prologue+num+epilogue, 0);
- n_mixed = vector<unsigned int>(prologue+num+epilogue, 0);
+ vector<unsigned int> n_static(prologue+num+epilogue, 0), n_forward(prologue+num+epilogue, 0),
+ n_backward(prologue+num+epilogue, 0), n_mixed(prologue+num+epilogue, 0);
for (int i = 0; i < static_cast<int>(prologue); i++)
if (variable_lag_lead[tmp_variable_reordered[i]].first != 0 && variable_lag_lead[tmp_variable_reordered[i]].second != 0)
@@ -1042,6 +1039,8 @@ ModelTree::computeBlockDecompositionAndFeedbackVariablesForEachBlock(const jacob
inv_variable_reordered[variable_reordered[i]] = i;
inv_equation_reordered[equation_reordered[i]] = i;
}
+
+ return { blocks, equation_lag_lead, variable_lag_lead, n_static, n_forward, n_backward, n_mixed };
}
void
@@ -1074,7 +1073,7 @@ ModelTree::printBlockDecomposition(const vector<pair<int, int>> &blocks) const
}
block_type_firstequation_size_mfs_t
-ModelTree::reduceBlocksAndTypeDetermination(const dynamic_jacob_map_t &dynamic_jacobian, vector<pair<int, int>> &blocks, const equation_type_and_normalized_equation_t &Equation_Type, const vector<int> &variable_reordered, const vector<int> &equation_reordered, vector<unsigned int> &n_static, vector<unsigned int> &n_forward, vector<unsigned int> &n_backward, vector<unsigned int> &n_mixed, vector<tuple<int, int, int, int>> &block_col_type, bool linear_decomposition)
+ModelTree::reduceBlocksAndTypeDetermination(const vector<pair<int, int>> &blocks, const equation_type_and_normalized_equation_t &Equation_Type, const vector<unsigned int> &n_static, const vector<unsigned int> &n_forward, const vector<unsigned int> &n_backward, const vector<unsigned int> &n_mixed, bool linear_decomposition)
{
int i = 0;
int count_equ = 0, blck_count_simult = 0;
@@ -1232,7 +1231,7 @@ ModelTree::reduceBlocksAndTypeDetermination(const dynamic_jacob_map_t &dynamic_j
}
vector<bool>
-ModelTree::equationLinear(map<tuple<int, int, int>, expr_t> first_order_endo_derivatives) const
+ModelTree::equationLinear(const map<tuple<int, int, int>, expr_t> &first_order_endo_derivatives) const
{
vector<bool> is_linear(symbol_table.endo_nbr(), true);
for (const auto &it : first_order_endo_derivatives)
@@ -1250,7 +1249,7 @@ ModelTree::equationLinear(map<tuple<int, int, int>, expr_t> first_order_endo_der
}
vector<bool>
-ModelTree::BlockLinear(const blocks_derivatives_t &blocks_derivatives, const vector<int> &variable_reordered) const
+ModelTree::BlockLinear() const
{
unsigned int nb_blocks = getNbBlocks();
vector<bool> blocks_linear(nb_blocks, true);
diff --git a/src/ModelTree.hh b/src/ModelTree.hh
index 90c896480083fd30fa44d1d4add6c2c1b7e46484..a9a32e99ec20c98220bed499ecb5a98eb82e281f 100644
--- a/src/ModelTree.hh
+++ b/src/ModelTree.hh
@@ -53,7 +53,7 @@ vectorToTuple(const vector<T> &v)
using equation_type_and_normalized_equation_t = vector<pair<EquationType, expr_t>>;
//! Vector describing variables: max_lag in the block, max_lead in the block
-using lag_lead_vector_t = vector<pair< int, int>>;
+using lag_lead_vector_t = vector<pair<int, int>>;
//! for each block contains tuple<Simulation_Type, first_equation, Block_Size, Recursive_part_Size>
using block_type_firstequation_size_mfs_t = vector<tuple<BlockSimulationType, int, int, int>>;
@@ -279,34 +279,38 @@ protected:
If no matching is found using a strictly positive cutoff, then a zero cutoff is applied (i.e. use a symbolic normalization); in that case, the method adds zeros in the jacobian matrices to reflect all the edges in the symbolic incidence matrix.
If no matching is found with a zero cutoff close to zero an error message is printout.
*/
- void computeNonSingularNormalization(jacob_map_t &contemporaneous_jacobian, double cutoff, jacob_map_t &static_jacobian, dynamic_jacob_map_t &dynamic_jacobian);
+ void computeNonSingularNormalization(jacob_map_t &contemporaneous_jacobian, double cutoff, jacob_map_t &static_jacobian);
//! Try to find a natural normalization if all equations are matched to an endogenous variable on the LHS
bool computeNaturalNormalization();
//! Try to normalized each unnormalized equation (matched endogenous variable only on the LHS)
multimap<int, int> computeNormalizedEquations() const;
//! Evaluate the jacobian and suppress all the elements below the cutoff
- void evaluateAndReduceJacobian(const eval_context_t &eval_context, jacob_map_t &contemporaneous_jacobian, jacob_map_t &static_jacobian, dynamic_jacob_map_t &dynamic_jacobian, double cutoff, bool verbose);
+ /*! Returns a pair (contemporaneous_jacobian, static_jacobian). Also fills
+ dynamic_jacobian. */
+ pair<jacob_map_t, jacob_map_t> evaluateAndReduceJacobian(const eval_context_t &eval_context, double cutoff, bool verbose);
//! Select and reorder the non linear equations of the model
- vector<pair<int, int>> select_non_linear_equations_and_variables(vector<bool> is_equation_linear, const dynamic_jacob_map_t &dynamic_jacobian, vector<int> &equation_reordered, vector<int> &variable_reordered,
- vector<int> &inv_equation_reordered, vector<int> &inv_variable_reordered,
- lag_lead_vector_t &equation_lag_lead, lag_lead_vector_t &variable_lag_lead,
- vector<unsigned int> &n_static, vector<unsigned int> &n_forward, vector<unsigned int> &n_backward, vector<unsigned int> &n_mixed);
+ /*! Returns a tuple (blocks, equation_lag_lead, variable_lag_lead, n_static,
+ n_forward, n_backward, n_mixed) */
+ tuple<vector<pair<int, int>>, lag_lead_vector_t, lag_lead_vector_t, vector<unsigned int>, vector<unsigned int>, vector<unsigned int>, vector<unsigned int>> select_non_linear_equations_and_variables(const vector<bool> &is_equation_linear);
//! Search the equations and variables belonging to the prologue and the epilogue of the model
- void computePrologueAndEpilogue(const jacob_map_t &static_jacobian, vector<int> &equation_reordered, vector<int> &variable_reordered);
+ void computePrologueAndEpilogue(const jacob_map_t &static_jacobian);
//! Determine the type of each equation of model and try to normalized the unnormalized equation using computeNormalizedEquations
- equation_type_and_normalized_equation_t equationTypeDetermination(const map<tuple<int, int, int>, expr_t> &first_order_endo_derivatives, const vector<int> &Index_Var_IM, const vector<int> &Index_Equ_IM, int mfs) const;
+ equation_type_and_normalized_equation_t equationTypeDetermination(const map<tuple<int, int, int>, expr_t> &first_order_endo_derivatives, int mfs) const;
//! Compute the block decomposition and for a non-recusive block find the minimum feedback set
- void computeBlockDecompositionAndFeedbackVariablesForEachBlock(const jacob_map_t &static_jacobian, const dynamic_jacob_map_t &dynamic_jacobian, vector<int> &equation_reordered, vector<int> &variable_reordered, vector<pair<int, int>> &blocks, const equation_type_and_normalized_equation_t &Equation_Type, bool verbose_, bool select_feedback_variable, int mfs, vector<int> &inv_equation_reordered, vector<int> &inv_variable_reordered, lag_lead_vector_t &equation_lag_lead, lag_lead_vector_t &variable_lag_lead_t, vector<unsigned int> &n_static, vector<unsigned int> &n_forward, vector<unsigned int> &n_backward, vector<unsigned int> &n_mixed) const;
+ /*! Returns a tuple (blocks, equation_lag_lead, variable_lag_lead, n_static,
+ n_forward, n_backward, n_mixed) */
+ tuple<vector<pair<int, int>>, lag_lead_vector_t, lag_lead_vector_t, vector<unsigned int>, vector<unsigned int>, vector<unsigned int>, vector<unsigned int>> computeBlockDecompositionAndFeedbackVariablesForEachBlock(const jacob_map_t &static_jacobian, const equation_type_and_normalized_equation_t &Equation_Type, bool verbose_, bool select_feedback_variable);
//! Reduce the number of block merging the same type equation in the prologue and the epilogue and determine the type of each block
- block_type_firstequation_size_mfs_t reduceBlocksAndTypeDetermination(const dynamic_jacob_map_t &dynamic_jacobian, vector<pair<int, int>> &blocks, const equation_type_and_normalized_equation_t &Equation_Type, const vector<int> &variable_reordered, const vector<int> &equation_reordered, vector<unsigned int> &n_static, vector<unsigned int> &n_forward, vector<unsigned int> &n_backward, vector<unsigned int> &n_mixed, vector<tuple<int, int, int, int>> &block_col_type, bool linear_decomposition);
+ block_type_firstequation_size_mfs_t reduceBlocksAndTypeDetermination(const vector<pair<int, int>> &blocks, const equation_type_and_normalized_equation_t &Equation_Type, const vector<unsigned int> &n_static, const vector<unsigned int> &n_forward, const vector<unsigned int> &n_backward, const vector<unsigned int> &n_mixed, bool linear_decomposition);
//! Determine the maximum number of lead and lag for the endogenous variable in a bloc
- void getVariableLeadLagByBlock(const dynamic_jacob_map_t &dynamic_jacobian, const vector<int> &components_set, int nb_blck_sim, lag_lead_vector_t &equation_lead_lag, lag_lead_vector_t &variable_lead_lag, const vector<int> &equation_reordered, const vector<int> &variable_reordered) const;
+ /*! Returns a pair { equation_lead,lag, variable_lead_lag } */
+ pair<lag_lead_vector_t, lag_lead_vector_t> getVariableLeadLagByBlock(const vector<int> &components_set, int nb_blck_sim) const;
//! For each equation determine if it is linear or not
- vector<bool> equationLinear(map<tuple<int, int, int>, expr_t> first_order_endo_derivatives) const;
+ vector<bool> equationLinear(const map<tuple<int, int, int>, expr_t> &first_order_endo_derivatives) const;
//! Print an abstract of the block structure of the model
void printBlockDecomposition(const vector<pair<int, int>> &blocks) const;
//! Determine for each block if it is linear or not
- vector<bool> BlockLinear(const blocks_derivatives_t &blocks_derivatives, const vector<int> &variable_reordered) const;
+ vector<bool> BlockLinear() const;
//! Remove equations specified by exclude_eqs
vector<int> includeExcludeEquations(set<pair<string, string>> &eqs, bool exclude_eqs,
vector<BinaryOpNode *> &equations, vector<int> &equations_lineno,
diff --git a/src/StaticModel.cc b/src/StaticModel.cc
index 6feb67291ca41cc25e8d49a2e6c9ac58f8a9ba12..eec56f50292d7f9ac2605840eb883dc7cbfc4d07 100644
--- a/src/StaticModel.cc
+++ b/src/StaticModel.cc
@@ -1144,35 +1144,27 @@ StaticModel::computingPass(int derivsOrder, int paramsDerivsOrder, const eval_co
if (block)
{
- jacob_map_t contemporaneous_jacobian, static_jacobian;
- vector<unsigned int> n_static, n_forward, n_backward, n_mixed;
+ auto [contemporaneous_jacobian, static_jacobian] = evaluateAndReduceJacobian(eval_context, cutoff, false);
- // for each block contains pair<Size, Feddback_variable>
- vector<pair<int, int>> blocks;
+ computeNonSingularNormalization(contemporaneous_jacobian, cutoff, static_jacobian);
- evaluateAndReduceJacobian(eval_context, contemporaneous_jacobian, static_jacobian, dynamic_jacobian, cutoff, false);
+ computePrologueAndEpilogue(static_jacobian);
- computeNonSingularNormalization(contemporaneous_jacobian, cutoff, static_jacobian, dynamic_jacobian);
+ auto first_order_endo_derivatives = collect_first_order_derivatives_endogenous();
- computePrologueAndEpilogue(static_jacobian, equation_reordered, variable_reordered);
-
- map<tuple<int, int, int>, expr_t> first_order_endo_derivatives = collect_first_order_derivatives_endogenous();
-
- equation_type_and_normalized_equation = equationTypeDetermination(first_order_endo_derivatives, variable_reordered, equation_reordered, mfs);
+ equation_type_and_normalized_equation = equationTypeDetermination(first_order_endo_derivatives, mfs);
cout << "Finding the optimal block decomposition of the model ..." << endl;
- lag_lead_vector_t equation_lag_lead, variable_lag_lead;
-
- computeBlockDecompositionAndFeedbackVariablesForEachBlock(static_jacobian, dynamic_jacobian, equation_reordered, variable_reordered, blocks, equation_type_and_normalized_equation, false, false, mfs, inv_equation_reordered, inv_variable_reordered, equation_lag_lead, variable_lag_lead, n_static, n_forward, n_backward, n_mixed);
+ auto [blocks, equation_lag_lead, variable_lag_lead, n_static, n_forward, n_backward, n_mixed] = computeBlockDecompositionAndFeedbackVariablesForEachBlock(static_jacobian, equation_type_and_normalized_equation, false, false);
- block_type_firstequation_size_mfs = reduceBlocksAndTypeDetermination(dynamic_jacobian, blocks, equation_type_and_normalized_equation, variable_reordered, equation_reordered, n_static, n_forward, n_backward, n_mixed, block_col_type, false);
+ block_type_firstequation_size_mfs = reduceBlocksAndTypeDetermination(blocks, equation_type_and_normalized_equation, n_static, n_forward, n_backward, n_mixed, false);
printBlockDecomposition(blocks);
- computeChainRuleJacobian(blocks_derivatives);
+ computeChainRuleJacobian();
- blocks_linear = BlockLinear(blocks_derivatives, variable_reordered);
+ blocks_linear = BlockLinear();
collect_block_first_order_derivatives();
@@ -2176,11 +2168,11 @@ StaticModel::get_Derivatives(int block)
}
void
-StaticModel::computeChainRuleJacobian(blocks_derivatives_t &blocks_derivatives)
+StaticModel::computeChainRuleJacobian()
{
map<int, expr_t> recursive_variables;
unsigned int nb_blocks = getNbBlocks();
- blocks_derivatives = blocks_derivatives_t(nb_blocks);
+ blocks_derivatives.resize(nb_blocks);
for (unsigned int block = 0; block < nb_blocks; block++)
{
block_derivatives_equation_variable_laglead_nodeid_t tmp_derivatives;
@@ -2191,7 +2183,6 @@ StaticModel::computeChainRuleJacobian(blocks_derivatives_t &blocks_derivatives)
int block_nb_recursives = block_size - block_nb_mfs;
if (simulation_type == SOLVE_TWO_BOUNDARIES_COMPLETE || simulation_type == SOLVE_TWO_BOUNDARIES_SIMPLE)
{
- blocks_derivatives.push_back(block_derivatives_equation_variable_laglead_nodeid_t(0));
for (int i = 0; i < block_nb_recursives; i++)
{
if (getBlockEquationType(block, i) == E_EVALUATE_S)
@@ -2220,7 +2211,6 @@ StaticModel::computeChainRuleJacobian(blocks_derivatives_t &blocks_derivatives)
}
else
{
- blocks_derivatives.push_back(block_derivatives_equation_variable_laglead_nodeid_t(0));
for (int i = 0; i < block_nb_recursives; i++)
{
if (getBlockEquationType(block, i) == E_EVALUATE_S)
diff --git a/src/StaticModel.hh b/src/StaticModel.hh
index b62d84bc2f8b6f53f773755bd98da58b8659d126..ed3c6227c800cc7aa379a83923dd0bb0ed40635d 100644
--- a/src/StaticModel.hh
+++ b/src/StaticModel.hh
@@ -90,7 +90,7 @@ private:
//! return a map on the block jacobian
map<tuple<int, int, int, int, int>, int> get_Derivatives(int block);
//! Computes chain rule derivatives of the Jacobian w.r. to endogenous variables
- void computeChainRuleJacobian(blocks_derivatives_t &blocks_derivatives);
+ void computeChainRuleJacobian();
//! Collect only the first derivatives
map<tuple<int, int, int>, expr_t> collect_first_order_derivatives_endogenous();