support for constant-pressure reactors

Cantera/clib/src/ctreactor.cpp

@@ -2,6 +2,7 @@
 // Cantera includes
 #include "zeroD/Reactor.h"
 #include "zeroD/FlowReactor.h"
+#include "zeroD/ConstPressureReactor.h"
 #include "zeroD/ReactorNet.h"
 #include "zeroD/Reservoir.h"
 #include "zeroD/Wall.h"
@@ -74,6 +75,8 @@ extern "C" {
             r = new Reactor();
         else if (type == FlowReactorType)
             r = new FlowReactor();
+        else if (type == ConstPressureReactorType)
+            r = new ConstPressureReactor();
         else if (type == ReservoirType)
             r = new Reservoir();
         else 

Cantera/cxx/include/zerodim.h

@@ -5,6 +5,9 @@
 #include "kernel/zeroD/Reservoir.h"
 #include "kernel/zeroD/Wall.h"
 #include "kernel/zeroD/flowControllers.h"
+#include "kernel/zeroD/FlowReactor.h"
+#include "kernel/zeroD/ConstPressureReactor.h"
+
 using namespace CanteraZeroD;
 
 #endif

Cantera/python/Cantera/Reactor.py

@@ -417,6 +417,38 @@ class FlowReactor(ReactorBase):
         _cantera.flowReactor_setMassFlowRate(self.__reactor_id, mdot)
         
         
+class ConstPressureReactor(ReactorBase):
+    """
+    """
+    def __init__(self, contents = None, name = '',
+                 volume = 1.0, energy = 'on',
+                 mdot = -1.0, 
+                 verbose = 0):
+        """
+        contents - Reactor contents. If not specified, the reactor is
+        initially empty. In this case, call method 'insert' to specify
+        the contents.
+
+        name - Used only to identify this reactor in output. If not
+        specified, defaults to 'Reactor_n', where n is an integer
+        assigned in the order Reactor objects are created.
+
+        volume - Initial reactor volume. Defaults to 1 m^3.
+
+        energy - Set to 'on' or 'off'. If set to 'off', the energy
+        equation is not solved, and the temperature is held at its
+        initial value. The default in 'on'.
+
+        verbose - if set to a non-zero value, additional diagnostic
+        information will be printed.
+        """
+        global _reactorcount
+        if name == '':
+            name = 'ConstPressureReactor_'+`_reactorcount`
+        _reactorcount += 1
+        ReactorBase.__init__(self, contents = contents, name = name,
+                             volume = volume, energy = energy,
+                             verbose = verbose, type = 4)
             
 
 class Reservoir(ReactorBase):

Cantera/src/CVodesIntegrator.cpp

@@ -11,8 +11,10 @@
 #include <iostream>
 using namespace std;
 
+#define OLD_SUNDIALS
 
 // sundials includes
+#ifdef OLD_SUNDIALS
 #include <sundialstypes.h>
 #include <sundialsmath.h>
 #include <cvodes.h>
@@ -21,6 +23,16 @@ using namespace std;
 #include <cvspgmr.h>
 #include <nvector.h>
 #include <nvector_serial.h>
+#else
+#include <sundials_types.h>
+#include <sundials_math.h>
+#include <cvodes.h>
+#include <cvodes_dense.h>
+#include <cvodes_diag.h>
+#include <cvodes_spgmr.h>
+//#include <nvector.h>
+#include <nvector_serial.h>
+#endif
 
 inline static N_Vector nv(void* x) {
     return reinterpret_cast<N_Vector>(x);
@@ -109,7 +121,11 @@ namespace Cantera {
         if (m_cvode_mem) {
             if (m_np > 0) 
                 CVodeSensFree(m_cvode_mem);
+#ifdef OLD_SUNDIALS
             CVodeFree(m_cvode_mem);
+#else
+            CVodeFree(&m_cvode_mem);
+#endif
         }
       if (m_y) N_VDestroy_Serial(nv(m_y));
       if (m_abstol) N_VDestroy_Serial(nv(m_abstol));

Cantera/src/CVodesIntegrator.h

@@ -26,6 +26,7 @@
 #include <nvector.h>
 #include <nvector_serial.h>
 
+
 namespace Cantera {
 
     class FuncData;

Cantera/src/FuncEval.h

@@ -23,7 +23,7 @@ namespace Cantera {
      *  Classes derived from FuncEval evaluate the right-hand-side function
      * \f$ \vec{F}(t,\vec{y})\f$ in 
      * \f[
-     *  \dot\vec{y} = \vec{F}(t,\vec{y}).
+     *  \dot{\vec{y}} = \vec{F}(t,\vec{y}).
      * \f] 
      *  @ingroup odeGroup 
      */
@@ -50,7 +50,9 @@ namespace Cantera {
          */
         virtual void getInitialConditions(double t0, size_t leny, double* y)=0;
 
-        /** Number of equations. */
+        /** 
+         * Number of equations. 
+         */
         virtual int neq()=0;
 
         /// Number of parameters.

Cantera/src/Makefile.in

@@ -216,7 +216,9 @@ clean:
 	(if test -d SunWS_cache ; then \
            $(RM) -rf SunWS_cache ; \
          fi )
+ifeq (@WITH_REACTORS@, 1)
 	cd zeroD; @MAKE@ clean
+endif
 	cd oneD;  @MAKE@ clean
 	cd converters; @MAKE@ clean
 	cd transport; @MAKE@ clean
@@ -226,7 +228,9 @@ depends: $(DEPENDS)
 	cat *.d > .depends
 	$(RM) $(DEPENDS) 
 	cd oneD; @MAKE@ depends
+ifeq (@WITH_REACTORS@, 1)
 	cd zeroD; @MAKE@ depends
+endif
 	cd converters; @MAKE@ depends
 	cd transport; @MAKE@ depends
 ifeq (@COMPILE_CATHERMO@, 1)

Cantera/src/Mu0Poly.h

@@ -1,7 +1,8 @@
 /**
  *  @file Mu0Poly.h
- *
- * $Author$
+ */
+
+/* $Author$
  * $Revision$
  * $Date$
  */

Cantera/src/NasaPoly1.h

@@ -1,3 +1,6 @@
+#ifndef CT_NASAPOLY1_H
+#define CT_NASAPOLY1_H
+
 /**
  *  @file NasaPoly1.h
  */
@@ -10,9 +13,6 @@
 // Copyright 2001  California Institute of Technology
 
 
-#ifndef CT_NASAPOLY1_H
-#define CT_NASAPOLY1_H
-
 #include "global.h"
 #include "SpeciesThermoInterpType.h"
 

Cantera/src/ShomatePoly.h

@@ -25,7 +25,7 @@ namespace Cantera {
      * \hat c_p(T) = A + B t + C t^2 + D t^3 + \frac{E}{t^2}
      * \f]
      * \f[
-     * \hat h^0(T)} = A t + \frac{B t^2}{2} + \frac{C t^3}{3} 
+     * \hat h^0(T) = A t + \frac{B t^2}{2} + \frac{C t^3}{3} 
      + \frac{D t^4}{4}  - \frac{E}{t}  + F.
      * \f]
      * \f[

Cantera/src/StoichManager.h

@@ -43,7 +43,7 @@ namespace Cantera {
      * \f[
      *  \dot C_k = \sum_k \nu^{(p)}_{k,i} R_i
      * \f]
-     * where \f$ \nu^{(p)_{k,i}}$ is the product-side stoichiometric
+     * where \f$ \nu^{(p)_{k,i}} \f$ is the product-side stoichiometric
      * coefficient of species \a k in reaction \a i.
      * This could be done be straightforward matrix multiplication, but would be inefficient, since most of the matrix elements of \f$ \nu^{(p)}_{k,i} \f$ are zero. We could do better by using sparse-matrix algorithms to compute this product. 
 

Cantera/src/importCTML.h

@@ -54,13 +54,7 @@ namespace Cantera {
      */
     void checkRxnElementBalance(Kinetics& kin, 
         const ReactionData &rdata, doublereal errorTolerance = 1.0e-3);
-    /**
-     * Extract the rate coefficient for a reaction from the xml node, kf.
-     * kf should point to a XML element named "rateCoeff".
-     * rdata is the partially filled ReactionData object for the reaction.
-     * This function will fill in more fields in the ReactionData object.
-     * 
-     */
+
     void getRateCoefficient(const XML_Node& kf, kinetics_t& kin, 
 			    ReactionData& rdata, int negA);
 

Cantera/src/transport/MixTransport.cpp

@@ -133,7 +133,7 @@ namespace Cantera {
      * and 
      * \f[
      * \Phi_{k,j} = \frac{\left[1 
-     * + \sqrt\left(\frac{\mu_k}{\mu_j}\sqrt{\frac{M_j}{M_k}\right)}\right]^2}
+     * + \sqrt\left(\frac{\mu_k}{\mu_j}\sqrt{\frac{M_j}{M_k}}\right)\right]^2}
      * {\sqrt{8}\sqrt{1 + M_k/M_j}}
      * \f] 
      * @see updateViscosity_T();