From 0b922739a627381737f82bb259390495282f6d19 Mon Sep 17 00:00:00 2001
From: Houtan Bastani <houtan@dynare.org>
Date: Thu, 12 May 2016 12:02:34 +0200
Subject: [PATCH] preprocessor: only compute first order derivatives w.r.t.
 parameters with identification. closes #1187

---
 ComputingTasks.cc |  2 +-
 DynamicModel.cc   |  6 ++--
 DynamicModel.hh   |  2 +-
 ModFile.cc        | 17 +++++++---
 ModelTree.cc      | 83 +++++++++++++++++++++++++----------------------
 ModelTree.hh      |  5 +--
 StaticModel.cc    |  6 ++--
 StaticModel.hh    |  2 +-
 8 files changed, 68 insertions(+), 55 deletions(-)

diff --git a/ComputingTasks.cc b/ComputingTasks.cc
index a9427bbb..a47e8497 100644
--- a/ComputingTasks.cc
+++ b/ComputingTasks.cc
@@ -1253,7 +1253,7 @@ PlannerObjectiveStatement::getPlannerObjective() const
 void
 PlannerObjectiveStatement::computingPass()
 {
-  model_tree->computingPass(eval_context_t(), false, true, true, false, false, false);
+  model_tree->computingPass(eval_context_t(), false, true, true, none, false, false);
 }
 
 void
diff --git a/DynamicModel.cc b/DynamicModel.cc
index bf165209..56dcc070 100644
--- a/DynamicModel.cc
+++ b/DynamicModel.cc
@@ -3113,7 +3113,7 @@ DynamicModel::runTrendTest(const eval_context_t &eval_context)
 }
 
 void
-DynamicModel::computingPass(bool jacobianExo, bool hessian, bool thirdDerivatives, bool paramsDerivatives,
+DynamicModel::computingPass(bool jacobianExo, bool hessian, bool thirdDerivatives, FileOutputType paramsDerivatives,
                             const eval_context_t &eval_context, bool no_tmp_terms, bool block, bool use_dll,
                             bool bytecode, bool compute_xrefs)
 {
@@ -3148,10 +3148,10 @@ DynamicModel::computingPass(bool jacobianExo, bool hessian, bool thirdDerivative
       computeHessian(vars);
     }
 
-  if (paramsDerivatives)
+  if (paramsDerivatives != none)
     {
       cout << " - derivatives of Jacobian/Hessian w.r. to parameters" << endl;
-      computeParamsDerivatives();
+      computeParamsDerivatives(paramsDerivatives);
 
       if (!no_tmp_terms)
         computeParamsDerivativesTemporaryTerms();
diff --git a/DynamicModel.hh b/DynamicModel.hh
index 13f18192..88f883cb 100644
--- a/DynamicModel.hh
+++ b/DynamicModel.hh
@@ -212,7 +212,7 @@ public:
     \param eval_context evaluation context for normalization
     \param no_tmp_terms if true, no temporary terms will be computed in the dynamic files
   */
-  void computingPass(bool jacobianExo, bool hessian, bool thirdDerivatives, bool paramsDerivatives,
+  void computingPass(bool jacobianExo, bool hessian, bool thirdDerivatives, FileOutputType paramsDerivatives,
                      const eval_context_t &eval_context, bool no_tmp_terms, bool block, bool use_dll, bool bytecode, bool compute_xrefs);
   //! Writes model initialization and lead/lag incidence matrix to output
   void writeOutput(ostream &output, const string &basename, bool block, bool byte_code, bool use_dll, int order, bool estimation_present, bool compute_xrefs, bool julia) const;
diff --git a/ModFile.cc b/ModFile.cc
index 1776fd9b..972a5985 100644
--- a/ModFile.cc
+++ b/ModFile.cc
@@ -489,8 +489,11 @@ ModFile::computingPass(bool no_tmp_terms, FileOutputType output, bool compute_xr
 
 	  const bool static_hessian = mod_file_struct.identification_present
 	    || mod_file_struct.estimation_analytic_derivation;
-	  const bool paramsDerivatives = mod_file_struct.identification_present
-	    || mod_file_struct.estimation_analytic_derivation;
+          FileOutputType paramsDerivatives = none;
+          if (mod_file_struct.identification_present)
+              paramsDerivatives = first;
+          if (mod_file_struct.estimation_analytic_derivation)
+              paramsDerivatives = third;
 	  static_model.computingPass(global_eval_context, no_tmp_terms, static_hessian,
 				     false, paramsDerivatives, block, byte_code);
 	}
@@ -502,7 +505,7 @@ ModFile::computingPass(bool no_tmp_terms, FileOutputType output, bool compute_xr
 	  || mod_file_struct.calib_smoother_present)
 	{
 	  if (mod_file_struct.perfect_foresight_solver_present)
-	    dynamic_model.computingPass(true, false, false, false, global_eval_context, no_tmp_terms, block, use_dll, byte_code, compute_xrefs);
+	    dynamic_model.computingPass(true, false, false, none, global_eval_context, no_tmp_terms, block, use_dll, byte_code, compute_xrefs);
 	      else
 		{
 		  if (mod_file_struct.stoch_simul_present
@@ -523,12 +526,16 @@ ModFile::computingPass(bool no_tmp_terms, FileOutputType output, bool compute_xr
 		  bool thirdDerivatives = mod_file_struct.order_option == 3 
 		    || mod_file_struct.estimation_analytic_derivation
 		    || output == third;
-		  bool paramsDerivatives = mod_file_struct.identification_present || mod_file_struct.estimation_analytic_derivation;
+                  FileOutputType paramsDerivatives = none;
+                  if (mod_file_struct.identification_present)
+                      paramsDerivatives = first;
+                  if (mod_file_struct.estimation_analytic_derivation)
+                      paramsDerivatives = third;
 		  dynamic_model.computingPass(true, hessian, thirdDerivatives, paramsDerivatives, global_eval_context, no_tmp_terms, block, use_dll, byte_code, compute_xrefs);
 		}
 	    }
 	  else // No computing task requested, compute derivatives up to 2nd order by default
-	    dynamic_model.computingPass(true, true, false, false, global_eval_context, no_tmp_terms, block, use_dll, byte_code, compute_xrefs);
+	    dynamic_model.computingPass(true, true, false, none, global_eval_context, no_tmp_terms, block, use_dll, byte_code, compute_xrefs);
     }
 
   for (vector<Statement *>::iterator it = statements.begin();
diff --git a/ModelTree.cc b/ModelTree.cc
index 703654fd..518bd3fd 100644
--- a/ModelTree.cc
+++ b/ModelTree.cc
@@ -1657,8 +1657,10 @@ ModelTree::sparseHelper(int order, ostream &output, int row_nb, int col_nb, Expr
 }
 
 void
-ModelTree::computeParamsDerivatives()
+ModelTree::computeParamsDerivatives(FileOutputType paramsDerivatives)
 {
+  if (!(paramsDerivatives == first || paramsDerivatives == second || paramsDerivatives == third))
+    return;
   set<int> deriv_id_set;
   addAllParamDerivId(deriv_id_set);
 
@@ -1675,18 +1677,19 @@ ModelTree::computeParamsDerivatives()
           residuals_params_derivatives[make_pair(eq, param)] = d1;
         }
 
-      for (first_derivatives_t::const_iterator it2 = residuals_params_derivatives.begin();
-           it2 != residuals_params_derivatives.end(); it2++)
-        {
-          int eq = it2->first.first;
-          int param1 = it2->first.second;
-          expr_t d1 = it2->second;
+      if (paramsDerivatives == second || paramsDerivatives == third)
+        for (first_derivatives_t::const_iterator it2 = residuals_params_derivatives.begin();
+             it2 != residuals_params_derivatives.end(); it2++)
+          {
+            int eq = it2->first.first;
+            int param1 = it2->first.second;
+            expr_t d1 = it2->second;
 
-          expr_t d2 = d1->getDerivative(param);
-          if (d2 == Zero)
-            continue;
-          residuals_params_second_derivatives[make_pair(eq, make_pair(param1, param))] = d2;
-        }
+            expr_t d2 = d1->getDerivative(param);
+            if (d2 == Zero)
+              continue;
+            residuals_params_second_derivatives[make_pair(eq, make_pair(param1, param))] = d2;
+          }
 
       for (first_derivatives_t::const_iterator it2 = first_derivatives.begin();
            it2 != first_derivatives.end(); it2++)
@@ -1701,33 +1704,35 @@ ModelTree::computeParamsDerivatives()
           jacobian_params_derivatives[make_pair(eq, make_pair(var, param))] = d2;
         }
 
-      for (second_derivatives_t::const_iterator it2 = jacobian_params_derivatives.begin();
-           it2 != jacobian_params_derivatives.end(); it2++)
-        {
-          int eq = it2->first.first;
-          int var = it2->first.second.first;
-          int param1 = it2->first.second.second;
-          expr_t d1 = it2->second;
-
-          expr_t d2 = d1->getDerivative(param);
-          if (d2 == Zero)
-            continue;
-          jacobian_params_second_derivatives[make_pair(eq, make_pair(var, make_pair(param1, param)))] = d2;
-        }
-
-      for (second_derivatives_t::const_iterator it2 = second_derivatives.begin();
-           it2 != second_derivatives.end(); it2++)
-        {
-          int eq = it2->first.first;
-          int var1 = it2->first.second.first;
-          int var2 = it2->first.second.second;
-          expr_t d1 = it2->second;
-
-          expr_t d2 = d1->getDerivative(param);
-          if (d2 == Zero)
-            continue;
-          hessian_params_derivatives[make_pair(eq, make_pair(var1, make_pair(var2, param)))] = d2;
-        }
+      if (paramsDerivatives == second || paramsDerivatives == third)
+        for (second_derivatives_t::const_iterator it2 = jacobian_params_derivatives.begin();
+             it2 != jacobian_params_derivatives.end(); it2++)
+          {
+            int eq = it2->first.first;
+            int var = it2->first.second.first;
+            int param1 = it2->first.second.second;
+            expr_t d1 = it2->second;
+
+            expr_t d2 = d1->getDerivative(param);
+            if (d2 == Zero)
+              continue;
+            jacobian_params_second_derivatives[make_pair(eq, make_pair(var, make_pair(param1, param)))] = d2;
+          }
+
+      if (paramsDerivatives == third)
+        for (second_derivatives_t::const_iterator it2 = second_derivatives.begin();
+             it2 != second_derivatives.end(); it2++)
+          {
+            int eq = it2->first.first;
+            int var1 = it2->first.second.first;
+            int var2 = it2->first.second.second;
+            expr_t d1 = it2->second;
+
+            expr_t d2 = d1->getDerivative(param);
+            if (d2 == Zero)
+              continue;
+            hessian_params_derivatives[make_pair(eq, make_pair(var1, make_pair(var2, param)))] = d2;
+          }
     }
 }
 
diff --git a/ModelTree.hh b/ModelTree.hh
index 259c3962..41f30d3b 100644
--- a/ModelTree.hh
+++ b/ModelTree.hh
@@ -1,5 +1,5 @@
 /*
- * Copyright (C) 2003-2015 Dynare Team
+ * Copyright (C) 2003-2016 Dynare Team
  *
  * This file is part of Dynare.
  *
@@ -29,6 +29,7 @@ using namespace std;
 #include <ostream>
 
 #include "DataTree.hh"
+#include "ExtendedPreprocessorTypes.hh"
 
 //! Vector describing equations: BlockSimulationType, if BlockSimulationType == EVALUATE_s then a expr_t on the new normalized equation
 typedef vector<pair<EquationType, expr_t > > equation_type_and_normalized_equation_t;
@@ -176,7 +177,7 @@ protected:
   /*! \param vars the derivation IDs w.r. to which derive the 2nd derivatives */
   void computeThirdDerivatives(const set<int> &vars);
   //! Computes derivatives of the Jacobian and Hessian w.r. to parameters
-  void computeParamsDerivatives();
+  void computeParamsDerivatives(FileOutputType paramsDerivatives);
   //! Write derivative of an equation w.r. to a variable
   void writeDerivative(ostream &output, int eq, int symb_id, int lag, ExprNodeOutputType output_type, const temporary_terms_t &temporary_terms) const;
   //! Computes temporary terms (for all equations and derivatives)
diff --git a/StaticModel.cc b/StaticModel.cc
index 8890cf1b..a3b42f56 100644
--- a/StaticModel.cc
+++ b/StaticModel.cc
@@ -1047,7 +1047,7 @@ StaticModel::collect_first_order_derivatives_endogenous()
 }
 
 void
-StaticModel::computingPass(const eval_context_t &eval_context, bool no_tmp_terms, bool hessian, bool thirdDerivatives, bool paramsDerivatives, bool block, bool bytecode)
+StaticModel::computingPass(const eval_context_t &eval_context, bool no_tmp_terms, bool hessian, bool thirdDerivatives, FileOutputType paramsDerivatives, bool block, bool bytecode)
 {
   initializeVariablesAndEquations();
 
@@ -1095,10 +1095,10 @@ StaticModel::computingPass(const eval_context_t &eval_context, bool no_tmp_terms
       computeThirdDerivatives(vars);
     }
 
-  if (paramsDerivatives)
+  if (paramsDerivatives != none)
     {
       cout << " - derivatives of Jacobian/Hessian w.r. to parameters" << endl;
-      computeParamsDerivatives();
+      computeParamsDerivatives(paramsDerivatives);
 
       if (!no_tmp_terms)
         computeParamsDerivativesTemporaryTerms();
diff --git a/StaticModel.hh b/StaticModel.hh
index 69b35d64..d8ec9790 100644
--- a/StaticModel.hh
+++ b/StaticModel.hh
@@ -164,7 +164,7 @@ public:
     \param hessian whether 2nd derivatives w.r. to exo, exo_det and endo should be computed
     \param paramsDerivatives whether 2nd derivatives w.r. to a pair (endo/exo/exo_det, parameter) should be computed
   */
-  void computingPass(const eval_context_t &eval_context, bool no_tmp_terms, bool hessian, bool thirdDerivatices, bool paramsDerivatives, bool block, bool bytecode);
+  void computingPass(const eval_context_t &eval_context, bool no_tmp_terms, bool hessian, bool thirdDerivatices, FileOutputType paramsDerivatives, bool block, bool bytecode);
 
   //! Adds informations for simulation in a binary file for a block decomposed model
   void Write_Inf_To_Bin_File_Block(const string &static_basename, const string &bin_basename, const int &num,
-- 
GitLab