-

Cantera/src/PureFluidPhase.h

@@ -0,0 +1,223 @@
+/**
+ *  @file PureFluid.h
+ *
+ * Declares class PureFluid
+ */
+
+// Copyright 2001  California Institute of Technology
+
+#ifndef CT_EOS_TPX_H
+#define CT_EOS_TPX_H
+
+#include "ThermoPhase.h"
+#include "../../ext/tpx/Sub.h"
+#include "../../ext/tpx/utils.h"
+
+namespace Cantera {
+    
+    /// Class for single-component fluids
+    class PureFluid  : public ThermoPhase {
+
+    public:
+
+        PureFluid() : ThermoPhase(), m_sub(0) {}
+
+        virtual ~PureFluid() { delete m_sub; }
+        
+        virtual void setParameters(int n, doublereal* c) {
+            if (n == 1) {
+                int subflag = int(c[0]);
+                if (m_sub) delete m_sub;
+                m_sub = tpx::GetSub(subflag);
+                if (m_sub == 0) {
+                    throw CanteraError("PureFluid::setParameters",
+                        "could not create new substance object.");
+                }
+                m_subflag = subflag;
+                m_mw = m_sub->MolWt();
+                double cp0_R, h0_RT, s0_R, T0, p;
+                T0 = 298.15;
+                if (T0 < m_sub->Tcrit()) {
+                    m_sub->Set(tpx::TX, T0, 1.0);
+                    p = 0.01*m_sub->P();
+                }
+                else {
+                    p = 0.001*m_sub->Pcrit();
+                }
+                m_sub->Set(tpx::TP, T0, p);
+                
+                m_spthermo->update_one(0, T0, &cp0_R, &h0_RT, &s0_R);
+                double s_R = s0_R - log(p/refPressure());
+                m_sub->setStdState(h0_RT*GasConstant*298.15/m_mw,
+                    s_R*GasConstant/m_mw, T0, p);
+            }
+        }
+        
+        virtual int eosType() const { return cPureFluid; }
+
+        /**
+         * Mixture molar enthalpy. Units: J/mol. 
+         */
+        virtual doublereal enthalpy_mole() const {
+            setTPXState();
+            doublereal h = m_sub->h() * m_mw;
+            check();
+            return h;            
+        }
+        
+        /**
+         * Mixture molar internal energy. Units: J/mol. 
+         */
+        virtual doublereal intEnergy_mole() const {
+            setTPXState();
+            doublereal u = m_sub->u() * m_mw;
+            check();
+            return u;            
+        }
+
+        /**
+         * Mixture molar entropy. Units: J/mol/K. 
+         */
+        virtual doublereal entropy_mole() const {
+            setTPXState();
+            doublereal s = m_sub->s() * m_mw;
+            check();
+            return s;            
+        }
+
+        /**
+         * Mixture molar Gibbs function. Units: J/mol. 
+         */
+        virtual doublereal gibbs_mole() const {
+            setTPXState();
+            doublereal g = m_sub->g() * m_mw;
+            check();
+            return g;            
+        }
+
+        /**
+         * Mixture molar heat capacity at constant pressure.
+         * Units: J/mol/K. 
+         */
+        virtual doublereal cp_mole() const {
+            setTPXState();
+            doublereal cp = m_sub->cp() * m_mw;
+            check();
+            return cp;            
+        }
+
+        /**
+         * Mixture molar heat capacity at constant volume.
+         * Units: J/mol/K. 
+         */
+        virtual doublereal cv_mole() const {
+            setTPXState();
+            doublereal cv = m_sub->cv() * m_mw;
+            check();
+            return cv;
+        }
+        
+        
+        /**
+         * Pressure. Units: Pa
+         */
+        virtual doublereal pressure() const {
+            setTPXState();
+            doublereal p = m_sub->P();
+            check();
+            return p;
+        }
+        
+        /**
+         * Set the pressure, holding temperature and composition
+         * fixed.  
+         */
+        virtual void setPressure(doublereal p) {
+            m_sub->Set(tpx::TP, temperature(), p);
+            setDensity(1.0/m_sub->v());
+            check();
+        }
+
+        virtual void getChemPotentials(doublereal* mu) const {
+            mu[0] = gibbs_mole();
+        }
+
+
+        tpx::Substance& TPX_Substance() { return *m_sub; }
+
+        /// critical temperature 
+        virtual doublereal critTemperature() const { return m_sub->Tcrit(); }
+        
+        /// critical pressure
+        virtual doublereal critPressure() const { return m_sub->Pcrit(); }
+        
+        /// critical density
+        virtual doublereal critDensity() const { return 1.0/m_sub->Vcrit(); }
+        
+        
+        /// saturation temperature
+        virtual doublereal satTemperature(doublereal p) const { 
+            doublereal ts = m_sub->Tsat(p);
+            check();
+            return ts;
+        }
+        
+        /// saturation pressure
+        virtual doublereal satPressure(doublereal t) const {
+            doublereal tsv = m_sub->Temp();
+            doublereal vsv = m_sub->v();
+            m_sub->Set(tpx::TP, t, 0.5*m_sub->Pcrit());
+            doublereal ps = m_sub->Ps();
+            m_sub->Set(tpx::TV,tsv,vsv);
+            check();
+            return ps;
+            }
+        
+        virtual doublereal vaporFraction() const {
+            setTPXState();
+            doublereal x = m_sub->x();
+            check();
+            return x;
+        }
+        
+        virtual void setState_satLiquid() {
+            setTPXState();
+            m_sub->Set(tpx::TX, temperature(), 0.0);
+            setDensity(1.0/m_sub->v());            
+            check();            
+        }
+        
+        virtual void setState_satVapor() {        
+            setTPXState();
+            m_sub->Set(tpx::TX, temperature(), 1.0);
+            setDensity(1.0/m_sub->v());
+            check();
+        }        
+
+protected:
+        
+        void setTPXState() const {
+            m_sub->Set(tpx::TV, temperature(), 1.0/density());            
+        }
+        
+        void check() const {
+            if (m_sub->Error()) {
+                throw CanteraError("PureFluidPhase",string(tpx::errorMsg(
+                                                               m_sub->Error())));
+            }
+        }
+        
+private:
+        mutable tpx::Substance* m_sub;
+        int m_subflag;
+        doublereal m_mw;
+    };
+
+}
+        
+#endif
+
+
+
+
+

Cantera/src/ThermoPhase.cpp

@@ -23,6 +23,67 @@
 
 namespace Cantera {
 
+    void ThermoPhase::setState_TPX(doublereal t, doublereal p, 
+        const doublereal* x) {
+        setMoleFractions(x); setTemperature(t); setPressure(p);
+    }
+
+    void ThermoPhase::setState_TPX(doublereal t, doublereal p, 
+        compositionMap& x) {
+        setMoleFractionsByName(x); setTemperature(t); setPressure(p);
+    }
+
+    void ThermoPhase::setState_TPX(doublereal t, doublereal p, 
+        const string& x) {
+        compositionMap xx;
+        int kk = nSpecies();
+        for (int k = 0; k < kk; k++) xx[speciesName(k)] = -1.0;
+        try {
+            parseCompString(x, xx);
+        }
+        catch (CanteraError) {
+            throw CanteraError("setState_TPX",
+                "Unknown species in composition map: "+ x);
+        }
+        setMoleFractionsByName(xx); setTemperature(t); setPressure(p);
+    }        
+
+    void ThermoPhase::setState_TPY(doublereal t, doublereal p, 
+        const doublereal* y) {
+        setMassFractions(y); setTemperature(t); setPressure(p);
+    }
+
+    void ThermoPhase::setState_TPY(doublereal t, doublereal p, 
+        compositionMap& y) {
+        setMassFractionsByName(y); setTemperature(t); setPressure(p);
+    }
+        
+    void ThermoPhase::setState_TPY(doublereal t, doublereal p, 
+        const string& y) {
+        compositionMap yy;
+        int kk = nSpecies();
+        for (int k = 0; k < kk; k++) yy[speciesName(k)] = -1.0;
+        try {
+            parseCompString(y, yy);
+        }
+        catch (CanteraError) {
+            throw CanteraError("setState_TPY",
+                "Unknown species in composition map: "+ y);
+        }
+        setMassFractionsByName(yy); setTemperature(t); setPressure(p);
+    }
+
+    void ThermoPhase::setState_TP(doublereal t, doublereal p) {
+        setTemperature(t); setPressure(p);
+    }
+
+    void ThermoPhase::setState_PX(doublereal p, doublereal* x) {
+        setMoleFractions(x); setPressure(p);
+    }
+
+    void ThermoPhase::setState_PY(doublereal p, doublereal* y) {
+        setMassFractions(y); setPressure(p);
+    }
 
     void ThermoPhase::setState_HP(doublereal h, doublereal p, 
         doublereal tol) {

Cantera/src/ThermoPhase.h

@@ -53,13 +53,14 @@ namespace Cantera {
     class ThermoPhase : public Phase {
 
     public:
+        
+        /// Constructor.
+        ThermoPhase() : Phase(), m_spthermo(0), m_speciesData(0),
+            m_index(-1), m_phi(0.0) {}
 
+        
         virtual ~ThermoPhase() {
             delete m_spthermo;
-	    // Taking this out because I think i don't own it
-            //delete m_speciesData;
-            m_spthermo = 0;
-            m_speciesData = 0;
         }
 
 
@@ -90,11 +91,19 @@ namespace Cantera {
         void setIndex(int m) { m_index = m; }
 
 
+        /// used to access data needed to construct transport manager
+        /// later.
         void saveSpeciesData(const XML_Node* data) {
             m_speciesData = data;
         }
 
-        const XML_Node* speciesData() { return m_speciesData; }
+        const XML_Node* speciesData() { 
+            if (m_speciesData) 
+                return m_speciesData;
+            else
+                throw CanteraError("ThermoPhase::speciesData",
+                    "m_speciesData is NULL");
+        }
 
 
         /**
@@ -236,11 +245,7 @@ namespace Cantera {
         }
 
         void setElectricPotential(doublereal v) {
-            //int nsp = nSpecies();
             m_phi = v;
-            //for (int k = 0; k < nsp; k++) {
-            //    setPotentialEnergy(k, v*charge(k)*Faraday);
-            //}
         }
 
         doublereal electricPotential() { return m_phi; }
@@ -484,63 +489,32 @@ namespace Cantera {
          * @{
          */
         /** Set the temperature (K), pressure (Pa), and mole fractions.  */
-        void setState_TPX(doublereal t, doublereal p, const doublereal* x) {
-            setMoleFractions(x); setTemperature(t); setPressure(p);
-        }
+        void setState_TPX(doublereal t, doublereal p, const doublereal* x);
 
         /** Set the temperature (K), pressure (Pa), and mole fractions.  */
-        void setState_TPX(doublereal t, doublereal p, compositionMap& x) {
-            setMoleFractionsByName(x); setTemperature(t); setPressure(p);
-        }
+        void setState_TPX(doublereal t, doublereal p, compositionMap& x);
 
         /** Set the temperature (K), pressure (Pa), and mole fractions.  */
-        void setState_TPX(doublereal t, doublereal p, const string& x) {
-            compositionMap xx;
-            int kk = nSpecies();
-            for (int k = 0; k < kk; k++) xx[speciesName(k)] = -1.0;
-	    try {
-	      parseCompString(x, xx);
-	    }
-            catch (CanteraError) {
-	      throw CanteraError("setState_TPX",
-				 "Unknown species in composition map: "+ x);
-	    }
-            setMoleFractionsByName(xx); setTemperature(t); setPressure(p);
-        }        
+        void setState_TPX(doublereal t, doublereal p, const string& x);
 
         /** Set the temperature (K), pressure (Pa), and mass fractions. */
-        void setState_TPY(doublereal t, doublereal p, const doublereal* y) {
-            setMassFractions(y); setTemperature(t); setPressure(p);
-        }
+        void setState_TPY(doublereal t, doublereal p, const doublereal* y);
 
         /** Set the temperature (K), pressure (Pa), and mass fractions. */
-        void setState_TPY(doublereal t, doublereal p, compositionMap& y) {
-            setMassFractionsByName(y); setTemperature(t); setPressure(p);
-        }
+        void setState_TPY(doublereal t, doublereal p, compositionMap& y);
         
         /** Set the temperature (K), pressure (Pa), and mass fractions.  */
-        void setState_TPY(doublereal t, doublereal p, const string& y) {
-            compositionMap yy;
-            int kk = nSpecies();
-            for (int k = 0; k < kk; k++) yy[speciesName(k)] = -1.0;
-            parseCompString(y, yy);
-            setMassFractionsByName(yy); setTemperature(t); setPressure(p);
-        }
+        void setState_TPY(doublereal t, doublereal p, const string& y);
 
         /** Set the temperature (K) and pressure (Pa) */
-        void setState_TP(doublereal t, doublereal p) {
-            setTemperature(t); setPressure(p);
-        }
+        void setState_TP(doublereal t, doublereal p);
 
         /** Set the pressure (Pa) and mole fractions.  */
-        void setState_PX(doublereal p, doublereal* x) {
-            setMoleFractions(x); setPressure(p);
-        }
+        void setState_PX(doublereal p, doublereal* x);
 
         /** Set the pressure (Pa) and mass fractions.  */
-        void setState_PY(doublereal p, doublereal* y) {
-            setMassFractions(y); setPressure(p);
-        }
+        void setState_PY(doublereal p, doublereal* y);
+
 
         /** Set the specific enthalpy (J/kg) and pressure (Pa). */
         void setState_HP(doublereal h, doublereal p, doublereal tol = 1.e-8);
@@ -577,7 +551,7 @@ namespace Cantera {
         void getActivities(doublereal* a) {
             getActivityConcentrations(a);
             int nsp = nSpecies();
-            doublereal rc = standardConcentration();
+            doublereal rc = 1.0/standardConcentration();
             scale(a, a + nsp, a, rc); 
         }
 
@@ -661,21 +635,14 @@ namespace Cantera {
         }
 
         doublereal minTemp(int k = -1) {
-            return m_spthermo->minTemp();
+            return m_spthermo->minTemp(k);
         }
 
         doublereal maxTemp(int k = -1) {
-            return m_spthermo->maxTemp();
+            return m_spthermo->maxTemp(k);
         }
         
             
-        ThermoPhase() {
-            m_spthermo = 0;
-	    m_speciesData = 0;
-            m_index = -1;
-            m_phi = 0.0;
-        }
-
     protected:
 
         /// Pointer to the species thermodynamic property manager

Cantera/src/transport/DustyGasTransport.cpp

@@ -7,7 +7,8 @@
  *
  */
 
-/*  $Author$
+/*
+ *  $Author$
  *  $Date$
  *  $Revision$
  *
@@ -24,24 +25,11 @@
 #endif
 
 #include "DustyGasTransport.h"
-#include "ctlapack.h"
-#include "../../ext/math/gmres.h"
-
-#include "DenseMatrix.h"
-#include "polyfit.h"
-#include "utilities.h"
-#include "TransportParams.h"
-#include "IdealGasPhase.h"
-
-#include "TransportFactory.h"
-
-#include <iostream>
 
 /** 
  * Mole fractions below MIN_X will be set to MIN_X when computing
  * transport properties.
  */
-
 #define MIN_X 1.e-20
 
 
@@ -115,28 +103,24 @@ namespace Cantera {
     void DustyGasTransport::updateBinaryDiffCoeffs() {
         if (m_bulk_ok) return;
         int n,m;
+
         // get the gaseous binary diffusion coefficients
-        //cout << "Gas binary diffusion coefficients: " << endl;
         m_gastran->getBinaryDiffCoeffs(m_nsp, m_d.begin());
         doublereal por2tort = m_porosity / m_tortuosity;
         for (n = 0; n < m_nsp; n++) 
             for (m = 0; m < m_nsp; m++) 
                 m_d(n,m) *= por2tort;
         m_bulk_ok = true;
-        //cout << m_d << endl;
     }
 
     void DustyGasTransport::updateKnudsenDiffCoeffs() {
         if (m_knudsen_ok) return;
         doublereal K_g = m_pore_radius * m_porosity / m_tortuosity;
         const doublereal FourThirds = 4.0/3.0;
-        //cout << "Knudsen diffusion coefficients: " << endl;
         for (int k = 0; k < m_nsp; k++) {
             m_dk[k] = FourThirds * K_g * sqrt((8.0 * GasConstant * m_temp)/
                 (Pi * m_mw[k]));
-            //cout << m_dk[k] << ", "; 
         }
-        //cout << endl;
         m_knudsen_ok = true;
     }
 
@@ -155,7 +139,6 @@ namespace Cantera {
             sum = 0.0;
             for (j = 0; j < m_nsp; j++) if (j != k) sum += m_x[j]/m_d(k,j);
             m_multidiff(k,k) = 1.0/m_dk[k] + sum;
-            //cout << "H matrix = " << endl << m_multidiff << endl;
         }
     }
 
@@ -197,8 +180,6 @@ namespace Cantera {
         // invert H
         int ierr = invert(m_multidiff);
 
-        //cout << "Diffusion coeff matrix: " << endl;
-        //cout << m_multidiff << endl;
         if (ierr != 0) {
             throw CanteraError("DustyGasTransport::updateMultiDiffCoeffs",
                 "invert returned ierr = "+int2str(ierr));

Cantera/src/transport/DustyGasTransport.h

@@ -1,9 +1,10 @@
-/**
- *
- *  @file DustyGasTransport.h
- *  Interface for class DustyGasTransport
- *
- */
+///
+///
+///  @file DustyGasTransport.h
+///  Interface for class DustyGasTransport
+///
+///
+
 
 // Copyright 2003  California Institute of Technology
 
@@ -11,23 +12,6 @@
 #ifndef CT_DUSTYGASTRAN_H
 #define CT_DUSTYGASTRAN_H
 
-
-// turn off warnings under Windows
-#ifdef WIN32
-#pragma warning(disable:4786)
-#pragma warning(disable:4503)
-#endif
-
-
-// STL includes
-#include <vector>
-#include <string>
-#include <map>
-#include <numeric>
-#include <algorithm>
-
-using namespace std;
-
 // Cantera includes
 #include "TransportBase.h"
 #include "../DenseMatrix.h"
@@ -35,45 +19,46 @@ using namespace std;
 
 namespace Cantera {
 
-
-    /**
-     * Class DustyGasTransport implements the Dusty Gas model for
-     * transport in porous media. As implemented here, only species
-     * transport is handled. The viscosity, thermal conductivity, and
-     * thermal diffusion coefficients are not implemented.
-    */
+    ///
+    /// Class DustyGasTransport implements the Dusty Gas model for
+    /// transport in porous media. As implemented here, only species
+    /// transport is handled. The viscosity, thermal conductivity, and
+    /// thermal diffusion coefficients are not implemented.
+    ///
     class DustyGasTransport : public Transport {
 
     public:
 
         /// default constructor
-      DustyGasTransport(thermo_t* thermo=0);
+        DustyGasTransport(thermo_t* thermo=0);
+        
+        /// Destructor. Does nothing, since class allocates no memory
+        /// on the heap.
         virtual ~DustyGasTransport() {}
-
+        
+        
+        //---------------------------------------------------------
         // overloaded base class methods
 
         virtual int model() { return cDustyGasTransport; }
 
         virtual void setParameters(int type, int k, doublereal* p);
-
-
-        //virtual void getBinaryDiffCoeffs(int ld, doublereal* d);
-
-        /**
-         * Get the multicomponent effective diffusion coefficients.
-         */
+        
         virtual void getMultiDiffCoeffs(int ld, doublereal* d);
-
-
+        
+        
+        
+        //-----------------------------------------------------------
         // new methods added in this class
-
-        /**
-         * Get the molar gas species fluxes. These fluxes include both the ordinary mass diffusion component
-         * and the Darcy (pressure-driven) commponent. 
-         */
+        
+        
+        /// Get the molar gas species fluxes. These fluxes include
+        /// both the ordinary mass diffusion component
+        /// and the Darcy (pressure-driven) commponent.
         void getMolarFluxes(const double* grad_conc,
             double grad_P, double* fluxes);
-
+        
+        
         /// Set the porosity (dimensionless)
         void setPorosity(doublereal porosity) {
             m_porosity = porosity;
@@ -98,27 +83,27 @@ namespace Cantera {
         void setMeanParticleDiameter(doublereal dbar) {
             m_diam = dbar;
         }
-
-        /** 
-         * Set the permeability. If not set, the value for
-         * close-packed spheres will be used by default. 
-         */ 
+        
+        /// Set the permeability. If not set, the value for
+        /// close-packed spheres will be used by default. 
         void setPermeability(doublereal B) {
             m_perm = B;
         }
-
+        
+        /// Return a reference to the transport manager used to compute the gas
+        /// binary diffusion coefficients and the visdcosity.
         Transport& gasTransport() { return *m_gastran; }
-        /**
-         * @internal
-         */
-
+        
+        
         friend class TransportFactory;
-
+        
+        
     protected:
 
         // called by TransportFactory
         void initialize(ThermoPhase* phase, Transport* gastr);
-
+        
+        
     private:
 
         void updateTransport_T();
@@ -147,10 +132,10 @@ namespace Cantera {
         vector_fp                    m_dk;
 
         /// temperature
-        doublereal m_temp;
+        doublereal                   m_temp;
 
         /// multicomponent diffusion coefficients
-        DenseMatrix  m_multidiff;
+        DenseMatrix                  m_multidiff;
 
         // work space
         vector_fp  m_spwork;
@@ -158,13 +143,13 @@ namespace Cantera {
         bool m_knudsen_ok;
         bool m_bulk_ok;
 
-        doublereal m_porosity;
-        doublereal m_tortuosity;
-        doublereal m_pore_radius;
-        doublereal m_diam;
-        doublereal m_perm;
+        doublereal m_porosity;      /// porosity
+        doublereal m_tortuosity;    /// tortuosity
+        doublereal m_pore_radius;   /// pore radius (m)
+        doublereal m_diam;          /// particle diameter (m)
+        doublereal m_perm;          /// permeability
 
-        Transport* m_gastran;
+        Transport* m_gastran;       /// pointer to gas transport manager
 
     };
 }