From 032bb9c36381aa2dace1541af1af1dc0f30f3e8a Mon Sep 17 00:00:00 2001
From: Ray Speth <speth@mit.edu>
Date: Tue, 19 Apr 2016 16:08:52 -0400
Subject: [PATCH] [Reactor] Move sensitivity parameter vector into FuncEval

This eliminates the need for class FuncData, and makes it possible to
introduce sensitivity parameters which have a nominal values other than 1.0.
---
 include/cantera/numerics/CVodesIntegrator.h |  3 ---
 include/cantera/numerics/FuncEval.h         |  7 +++++-
 src/numerics/CVodesIntegrator.cpp           | 27 ++++-----------------
 src/zeroD/ReactorNet.cpp                    |  1 +
 4 files changed, 12 insertions(+), 26 deletions(-)

diff --git a/include/cantera/numerics/CVodesIntegrator.h b/include/cantera/numerics/CVodesIntegrator.h
index fed74b377..307adbdd5 100644
--- a/include/cantera/numerics/CVodesIntegrator.h
+++ b/include/cantera/numerics/CVodesIntegrator.h
@@ -14,8 +14,6 @@
 namespace Cantera
 {
 
-class FuncData;
-
 /**
  * Exception thrown when a CVODES error is encountered.
  * @deprecated Unused. To be removed after Cantera 2.3.
@@ -109,7 +107,6 @@ private:
     double m_hmax, m_hmin;
     int m_maxsteps;
     int m_maxErrTestFails;
-    std::unique_ptr<FuncData> m_fdata;
     N_Vector* m_yS;
     size_t m_np;
     int m_mupper, m_mlower;
diff --git a/include/cantera/numerics/FuncEval.h b/include/cantera/numerics/FuncEval.h
index bdef34159..a63bbc79e 100644
--- a/include/cantera/numerics/FuncEval.h
+++ b/include/cantera/numerics/FuncEval.h
@@ -58,8 +58,13 @@ public:
 
     //! Number of sensitivity parameters.
     virtual size_t nparams() {
-        return 0;
+        return m_sens_params.size();
     }
+
+    //! Values for the problem parameters for which sensitivities are computed
+    //! This is the array which is perturbed and passed back as the fourth
+    //! argument to eval().
+    vector_fp m_sens_params;
 };
 
 }
diff --git a/src/numerics/CVodesIntegrator.cpp b/src/numerics/CVodesIntegrator.cpp
index 31e81fb03..c118b56a7 100644
--- a/src/numerics/CVodesIntegrator.cpp
+++ b/src/numerics/CVodesIntegrator.cpp
@@ -33,18 +33,6 @@ typedef long int sd_size_t;
 namespace Cantera
 {
 
-class FuncData
-{
-public:
-    FuncData(FuncEval* f, size_t npar = 0) {
-        m_pars.resize(npar, 1.0);
-        m_func = f;
-    }
-    virtual ~FuncData() {}
-    vector_fp m_pars;
-    FuncEval* m_func;
-};
-
 extern "C" {
     /**
      * Function called by cvodes to evaluate ydot given y.  The CVODE integrator
@@ -54,13 +42,11 @@ extern "C" {
      * evaluates the desired equations.
      * @ingroup odeGroup
      */
-    static int cvodes_rhs(realtype t, N_Vector y, N_Vector ydot,
-                          void* f_data)
+    static int cvodes_rhs(realtype t, N_Vector y, N_Vector ydot, void* f_data)
     {
         try {
-            FuncData* d = (FuncData*)f_data;
-            FuncEval* f = d->m_func;
-            f->eval(t, NV_DATA_S(y), NV_DATA_S(ydot), d->m_pars.data());
+            FuncEval* f = (FuncEval*) f_data;
+            f->eval(t, NV_DATA_S(y), NV_DATA_S(ydot), f->m_sens_params.data());
         } catch (CanteraError& err) {
             std::cerr << err.what() << std::endl;
             return 1; // possibly recoverable error
@@ -309,17 +295,14 @@ void CVodesIntegrator::initialize(double t0, FuncEval& func)
         }
     }
 
-    // pass a pointer to func in m_data
-    m_fdata.reset(new FuncData(&func, func.nparams()));
-
-    flag = CVodeSetUserData(m_cvode_mem, m_fdata.get());
+    flag = CVodeSetUserData(m_cvode_mem, &func);
     if (flag != CV_SUCCESS) {
         throw CanteraError("CVodesIntegrator::initialize",
                            "CVodeSetUserData failed.");
     }
     if (func.nparams() > 0) {
         sensInit(t0, func);
-        flag = CVodeSetSensParams(m_cvode_mem, m_fdata->m_pars.data(),
+        flag = CVodeSetSensParams(m_cvode_mem, func.m_sens_params.data(),
                                   NULL, NULL);
     }
     applyOptions();
diff --git a/src/zeroD/ReactorNet.cpp b/src/zeroD/ReactorNet.cpp
index dccad6ef9..a0aded592 100644
--- a/src/zeroD/ReactorNet.cpp
+++ b/src/zeroD/ReactorNet.cpp
@@ -266,6 +266,7 @@ void ReactorNet::registerSensitivityReaction(void* reactor,
     m_paramNames.push_back(name);
     m_sensOrder[R][reactionIndex] = m_ntotpar;
     m_ntotpar++;
+    m_sens_params.push_back(1.0);
 }
 
 }