From d246f9f99a1c78fa7a76fe407a200103fed6f2ba Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?S=C3=A9bastien=20Villemot?= <sebastien@dynare.org>
Date: Mon, 24 Apr 2023 17:47:29 +0200
Subject: [PATCH] ModelTree::computeNormalization(): throw an exception in case
of normalization failure
It would previously return a boolean. The exception is more convenient for
producing a different error message in the case of the specialized algorithm
introduced in the next commit.
---
src/ModelTree.cc | 70 ++++++++++++++++++++++++------------------------
src/ModelTree.hh | 7 ++++-
2 files changed, 41 insertions(+), 36 deletions(-)
diff --git a/src/ModelTree.cc b/src/ModelTree.cc
index 37f844db..7eef0981 100644
--- a/src/ModelTree.cc
+++ b/src/ModelTree.cc
@@ -234,8 +234,8 @@ ModelTree::operator=(const ModelTree &m)
return *this;
}
-bool
-ModelTree::computeNormalization(const jacob_map_t &contemporaneous_jacobian, bool verbose)
+void
+ModelTree::computeNormalization(const jacob_map_t &contemporaneous_jacobian)
{
const int n {static_cast<int>(equations.size())};
@@ -261,19 +261,11 @@ ModelTree::computeNormalization(const jacob_map_t &contemporaneous_jacobian, boo
// Check if all variables are normalized
if (auto it = find(mate_map.begin(), mate_map.begin() + n, boost::graph_traits<BipartiteGraph>::null_vertex());
it != mate_map.begin() + n)
- {
- if (verbose)
- cerr << "Could not normalize the " << modelClassName() << ". Variable "
- << symbol_table.getName(symbol_table.getID(SymbolType::endogenous, it - mate_map.begin()))
- << " is not in the maximum cardinality matching." << endl;
- return false;
- }
+ throw ModelNormalizationFailed {symbol_table.getName(symbol_table.getID(SymbolType::endogenous, it - mate_map.begin())) };
// Create the resulting map, by copying the n first elements of mate_map, and substracting n to them
endo2eq.resize(equations.size());
transform(mate_map.begin(), mate_map.begin() + n, endo2eq.begin(), [=](vertex_descriptor_t i) { return i-n; });
-
- return true;
}
bool
@@ -307,9 +299,8 @@ ModelTree::computeNonSingularNormalization(const eval_context_t &eval_context)
double current_cutoff = 0.99999999;
const double cutoff_lower_limit = 1e-19;
- bool found_normalization = false;
int last_suppressed = 0;
- while (!found_normalization && current_cutoff > cutoff_lower_limit)
+ while (current_cutoff > cutoff_lower_limit)
{
// Drop elements below cutoff from normalized contemporaneous Jacobian
jacob_map_t normalized_contemporaneous_jacobian_above_cutoff;
@@ -321,35 +312,44 @@ ModelTree::computeNonSingularNormalization(const eval_context_t &eval_context)
suppressed++;
if (suppressed != last_suppressed)
- found_normalization = computeNormalization(normalized_contemporaneous_jacobian_above_cutoff, false);
+ try
+ {
+ computeNormalization(normalized_contemporaneous_jacobian_above_cutoff);
+ return true;
+ }
+ catch (ModelNormalizationFailed &e)
+ {
+ }
last_suppressed = suppressed;
- if (!found_normalization)
- {
- current_cutoff /= 2;
- // In this last case try to normalize with the complete jacobian
- if (current_cutoff <= cutoff_lower_limit)
- found_normalization = computeNormalization(normalized_contemporaneous_jacobian, false);
- }
+ current_cutoff /= 2;
}
- if (!found_normalization)
+ // Try to normalize with the complete numerical jacobian
+ try
+ {
+ computeNormalization(normalized_contemporaneous_jacobian);
+ return true;
+ }
+ catch (ModelNormalizationFailed &e)
{
- cout << "Normalization failed with cutoff, trying symbolic normalization..." << endl;
- /* If no non-singular normalization can be found, try to find a
- normalization even with a potential singularity. */
- auto symbolic_jacobian = computeSymbolicJacobian(true);
- found_normalization = computeNormalization(symbolic_jacobian, true);
}
-#ifdef DEBUG
- for (size_t i {0}; i < equations.size(); i++)
- cout << "Variable " << symbol_table.getName(symbol_table.getID(SymbolType::endogenous, i))
- << " associated to equation " << endo2eq[i] << " (" << equation_tags.getTagValueByEqnAndKey(endo2eq[i], "name") << ")" << endl;
-#endif
+ cout << "Normalization failed with cutoff, trying symbolic normalization..." << endl;
+ /* If no non-singular normalization can be found, try to find a
+ normalization even with a potential singularity. */
+ auto symbolic_jacobian { computeSymbolicJacobian(true) };
+ try
+ {
+ computeNormalization(symbolic_jacobian);
+ return true;
+ }
+ catch (ModelNormalizationFailed &e)
+ {
+ cerr << "Could not normalize the " << modelClassName() << ". Variable "
+ << e.unmatched_endo << " is not in the maximum cardinality matching." << endl;
+ }
- /* NB: If normalization failed, an explanatory message has been printed by the last call
- to computeNormalization(), which has verbose=true */
- return found_normalization;
+ return false;
}
ModelTree::jacob_map_t
diff --git a/src/ModelTree.hh b/src/ModelTree.hh
index 8d3aff50..7ba465a5 100644
--- a/src/ModelTree.hh
+++ b/src/ModelTree.hh
@@ -436,12 +436,17 @@ private:
// Compute {var,eq}_idx_orig2block from {var,eq}_idx_block2orig
void updateReverseVariableEquationOrderings();
+ struct ModelNormalizationFailed
+ {
+ string unmatched_endo; // Name of endogenous not in maximum cardinality matching
+ };
+
//! Compute the matching between endogenous and variable using the jacobian contemporaneous_jacobian
/*!
\param contemporaneous_jacobian Jacobian used as an incidence matrix: all elements declared in the map (even if they are zero), are used as vertices of the incidence matrix
\return True if a complete normalization has been achieved
*/
- bool computeNormalization(const jacob_map_t &contemporaneous_jacobian, bool verbose);
+ void computeNormalization(const jacob_map_t &contemporaneous_jacobian);
//! Try to compute the matching between endogenous and variable using a decreasing cutoff
/*!
--
GitLab