diff --git a/Cantera/src/base/xml.cpp b/Cantera/src/base/xml.cpp index 02b83c335..5c6396fbc 100755 --- a/Cantera/src/base/xml.cpp +++ b/Cantera/src/base/xml.cpp @@ -482,6 +482,27 @@ namespace Cantera { } } + void XML_Node::clear() { + int n = static_cast(m_children.size()); + for (int i = 0; i < n; i++) { + if (m_children[i]) { + if (m_children[i]->parent() == this) { + delete m_children[i]; + m_children[i] = 0; + } + } + } + m_value.clear(); + m_childindex.clear(); + m_attribs.clear(); + m_children.clear(); + + m_nchildren = 0; + m_iscomment = false; + m_linenum = 0; + + } + // Add a child node to the current node containing a comment /* * Child node will have the name, "comment". diff --git a/Cantera/src/base/xml.h b/Cantera/src/base/xml.h index 4f48c5505..db4d6ed51 100755 --- a/Cantera/src/base/xml.h +++ b/Cantera/src/base/xml.h @@ -380,6 +380,13 @@ namespace Cantera { */ std::string attrib(const std::string & attr) const; + //! Clear the current node and everything under it + /*! + * The value, attributes and children are all zeroed. The name and the + * parent information is kept. + */ + void clear(); + private: //! Returns a changeable value of the attributes map for the current node /*! @@ -644,7 +651,7 @@ namespace Cantera { //! into the destination XML_Node tree, doing a complete copy //! as we go. /*! - * Note this is a const function becuase the current XML_Node and + * Note this is a const function because the current XML_Node and * its children isn't altered by this operation. * * @param node_dest This is the XML node to receive the information diff --git a/Cantera/src/kinetics/Enhanced3BConc.h b/Cantera/src/kinetics/Enhanced3BConc.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/FalloffFactory.cpp b/Cantera/src/kinetics/FalloffFactory.cpp old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/FalloffFactory.h b/Cantera/src/kinetics/FalloffFactory.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/FalloffMgr.h b/Cantera/src/kinetics/FalloffMgr.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/GRI_30_Kinetics.cpp b/Cantera/src/kinetics/GRI_30_Kinetics.cpp old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/GRI_30_Kinetics.h b/Cantera/src/kinetics/GRI_30_Kinetics.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/GasKinetics.cpp b/Cantera/src/kinetics/GasKinetics.cpp old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/GasKinetics.h b/Cantera/src/kinetics/GasKinetics.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/GasKineticsWriter.cpp b/Cantera/src/kinetics/GasKineticsWriter.cpp old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/GasKineticsWriter.h b/Cantera/src/kinetics/GasKineticsWriter.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/Group.cpp b/Cantera/src/kinetics/Group.cpp old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/Group.h b/Cantera/src/kinetics/Group.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/ImplicitChem.cpp b/Cantera/src/kinetics/ImplicitChem.cpp old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/ImplicitChem.h b/Cantera/src/kinetics/ImplicitChem.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/ImplicitSurfChem.cpp b/Cantera/src/kinetics/ImplicitSurfChem.cpp old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/ImplicitSurfChem.h b/Cantera/src/kinetics/ImplicitSurfChem.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/Kinetics.h b/Cantera/src/kinetics/Kinetics.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/RateCoeffMgr.h b/Cantera/src/kinetics/RateCoeffMgr.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/ReactionData.h b/Cantera/src/kinetics/ReactionData.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/ReactionPath.cpp b/Cantera/src/kinetics/ReactionPath.cpp old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/ReactionPath.h b/Cantera/src/kinetics/ReactionPath.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/RxnRates.h b/Cantera/src/kinetics/RxnRates.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/StoichManager.h b/Cantera/src/kinetics/StoichManager.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/ThirdBodyMgr.h b/Cantera/src/kinetics/ThirdBodyMgr.h old mode 100755 new mode 100644 diff --git a/Cantera/src/kinetics/reaction_defs.h b/Cantera/src/kinetics/reaction_defs.h old mode 100755 new mode 100644 diff --git a/Cantera/src/numerics/Integrator.h b/Cantera/src/numerics/Integrator.h old mode 100755 new mode 100644 index 47fc811d8..2006f62e8 --- a/Cantera/src/numerics/Integrator.h +++ b/Cantera/src/numerics/Integrator.h @@ -64,7 +64,14 @@ namespace Cantera { public: - virtual ~Integrator() {} + //! Default Constructor + Integrator() + { + } + + //! Destructor + virtual ~Integrator() { + } /** Set or reset the number of equations. */ //virtual void resize(int n)=0; @@ -76,10 +83,12 @@ namespace Cantera { * @param abstol array of N absolute tolerance values */ virtual void setTolerances(doublereal reltol, int n, - doublereal* abstol) { warn("setTolerances"); } + doublereal* abstol) { + warn("setTolerances"); + } - /** - * Set error tolerances. + //! Set error tolerances. + /*! * @param reltol scalar relative tolerance * @param abstol scalar absolute tolerance */ diff --git a/Cantera/src/thermo/ThermoFactory.cpp b/Cantera/src/thermo/ThermoFactory.cpp index 85ae9aacf..70ad4164d 100644 --- a/Cantera/src/thermo/ThermoFactory.cpp +++ b/Cantera/src/thermo/ThermoFactory.cpp @@ -418,6 +418,17 @@ namespace Cantera { ", is not a phase element."); } + /* + * In this section of code, we get the reference to the + * phase xml tree within the ThermoPhase object. Then, + * we clear it and fill it with the current information that + * we are about to use to construct the object. We will then + * be able to resurrect the information later by calling xml(). + */ + XML_Node &phaseNode_XML = th->xml(); + phaseNode_XML.clear(); + phase.copy(&phaseNode_XML); + // set the id attribute of the phase to the 'id' attribute // in the XML tree. th->setID(phase.id()); diff --git a/Cantera/src/transport/AqueousTransport.cpp b/Cantera/src/transport/AqueousTransport.cpp index abb947491..55e820462 100644 --- a/Cantera/src/transport/AqueousTransport.cpp +++ b/Cantera/src/transport/AqueousTransport.cpp @@ -19,6 +19,7 @@ #include "utilities.h" #include "TransportParams.h" +#include "LiquidTransportParams.h" #include "TransportFactory.h" #include "ctlapack.h" @@ -72,7 +73,7 @@ namespace Cantera { /* * This is where we dimension everything. */ - bool AqueousTransport::init(TransportParams& tr) { + bool AqueousTransport::initLiquid( LiquidTransportParams& tr ) { // constant substance attributes m_thermo = tr.thermo; @@ -86,15 +87,11 @@ namespace Cantera { m_thermo->molecularWeights().end(), m_mw.begin()); // copy polynomials and parameters into local storage - m_poly = tr.poly; m_visccoeffs = tr.visccoeffs; m_condcoeffs = tr.condcoeffs; m_diffcoeffs = tr.diffcoeffs; - m_mode = tr.mode; - m_diam = tr.diam; - m_eps = tr.eps; - m_alpha = tr.alpha; + m_mode = tr.mode_; m_phi.resize(m_nsp, m_nsp, 0.0); @@ -184,6 +181,7 @@ namespace Cantera { /******************* binary diffusion coefficients **************/ + //================================================================================================ void AqueousTransport::getBinaryDiffCoeffs(const int ld, doublereal* const d) { int i,j; @@ -200,32 +198,52 @@ namespace Cantera { d[ld*j + i] = rp * m_bdiff(i,j); } } - - + //================================================================================================ + // Get the electrical Mobilities (m^2/V/s). + /* + * This function returns the mobilities. In some formulations + * this is equal to the normal mobility multiplied by faraday's constant. + * + * Frequently, but not always, the mobility is calculated from the + * diffusion coefficient using the Einstein relation + * + * \f[ + * \mu^e_k = \frac{F D_k}{R T} + * \f] + * + * @param mobil_e Returns the mobilities of + * the species in array \c mobil_e. The array must be + * dimensioned at least as large as the number of species. + */ void AqueousTransport::getMobilities(doublereal* const mobil) { - // this needs to be checked out. int k; getMixDiffCoeffs(DATA_PTR(m_spwork)); doublereal c1 = ElectronCharge / (Boltzmann * m_temp); for (k = 0; k < m_nsp; k++) { - mobil[k] = c1 * m_spwork[k] * m_thermo->charge(k); + mobil[k] = c1 * m_spwork[k]; } } - - - + //================================================================================================ + void AqueousTransport::getFluidMobilities(doublereal* const mobil) { + getMixDiffCoeffs(DATA_PTR(m_spwork)); + doublereal c1 = 1.0 / (GasConstant * m_temp); + for (int k = 0; k < m_nsp; k++) { + mobil[k] = c1 * m_spwork[k]; + } + } + //================================================================================================ void AqueousTransport::set_Grad_V(const doublereal* const grad_V) { for (int a = 0; a < m_nDim; a++) { m_Grad_V[a] = grad_V[a]; } } - + //================================================================================================ void AqueousTransport::set_Grad_T(const doublereal* const grad_T) { for (int a = 0; a < m_nDim; a++) { m_Grad_T[a] = grad_T[a]; } } - + //================================================================================================ void AqueousTransport::set_Grad_X(const doublereal* const grad_X) { int itop = m_nDim * m_nsp; for (int i = 0; i < itop; i++) { @@ -585,17 +603,15 @@ namespace Cantera { * This function returns a Transport data object for a given species. * */ - struct GasTransportData AqueousTransport:: - getGasTransportData(int kSpecies) + struct LiquidTransportData AqueousTransport:: + getLiquidTransportData(int kSpecies) { - struct GasTransportData td; + struct LiquidTransportData td; td.speciesName = m_thermo->speciesName(kSpecies); - - td.wellDepth = m_eps[kSpecies] / Boltzmann; - td.diameter = m_diam(kSpecies, kSpecies) * 1.0E10; - td.polarizability = m_alpha[kSpecies] * 1.0E30; - + /* NEEDS WORK + td.hydroradius = ???; + */ return td; } diff --git a/Cantera/src/transport/AqueousTransport.h b/Cantera/src/transport/AqueousTransport.h index 8a460a2c9..1d593bd16 100644 --- a/Cantera/src/transport/AqueousTransport.h +++ b/Cantera/src/transport/AqueousTransport.h @@ -11,13 +11,15 @@ #ifndef CT_AQUEOUSTRAN_H -#define CT_AQYEOUSTRAN_H +#define CT_AQUEOUSTRAN_H using namespace std; // Cantera includes #include "TransportBase.h" #include "DenseMatrix.h" +#include "TransportParams.h" +#include "LiquidTransportParams.h" #include @@ -29,7 +31,7 @@ using namespace std; namespace Cantera { - class TransportParams; + class LiquidTransportParams; //! Class AqueousTransport implements mixture-averaged transport @@ -133,7 +135,7 @@ namespace Cantera { virtual ~AqueousTransport() {} //! Return the model id for this transport parameterization - virtual int model() { return cAqueousTransport; } + virtual int model() const { return cAqueousTransport; } //! overloaded base class methods @@ -192,12 +194,43 @@ namespace Cantera { */ virtual void getMixDiffCoeffs(doublereal* const d); - - //! Get the Mobilities + //! Get the Electrical mobilities (m^2/V/s). /*! - * @param mobil + * This function returns the electrical mobilities. In some formulations + * this is equal to the normal mobility multiplied by faraday's constant. + * + * Frequently, but not always, the mobility is calculated from the + * diffusion coefficient using the Einstein relation + * + * \f[ + * \mu^e_k = \frac{F D_k}{R T} + * \f] + * + * @param mobil_e Returns the mobilities of + * the species in array \c mobil_e. The array must be + * dimensioned at least as large as the number of species. */ - virtual void getMobilities(doublereal* const mobil); + virtual void getMobilities(doublereal* const mobil_e); + + //! Get the fluid mobilities (s kmol/kg). + /*! + * This function returns the fluid mobilities. Usually, you have + * to multiply Faraday's constant into the resulting expression + * to general a species flux expression. + * + * Frequently, but not always, the mobility is calculated from the + * diffusion coefficient using the Einstein relation + * + * \f[ + * \mu^f_k = \frac{D_k}{R T} + * \f] + * + * @param mobil_f Returns the mobilities of + * the species in array \c mobil_f. The array must be + * dimensioned at least as large as the number of species. + */ + virtual void getFluidMobilities(doublereal* const mobil_f); + //! Specify the value of the gradient of the voltage /*! @@ -275,7 +308,7 @@ namespace Cantera { * @param tr Transport parameters for all of the species * in the phase. */ - virtual bool init(TransportParams& tr); + virtual bool initLiquid( LiquidTransportParams& tr ); friend class TransportFactory; @@ -286,7 +319,7 @@ namespace Cantera { * * @param k Species number to obtain the properties about. */ - struct GasTransportData getGasTransportData(int k); + struct LiquidTransportData getLiquidTransportData(int k); //! Solve the stefan_maxell equations for the diffusive fluxes. diff --git a/Cantera/src/transport/DustyGasTransport.h b/Cantera/src/transport/DustyGasTransport.h index aaedc9947..e28fa4ae5 100644 --- a/Cantera/src/transport/DustyGasTransport.h +++ b/Cantera/src/transport/DustyGasTransport.h @@ -40,7 +40,7 @@ namespace Cantera { //--------------------------------------------------------- // overloaded base class methods - virtual int model() { return cDustyGasTransport; } + virtual int model() const { return cDustyGasTransport; } virtual void setParameters(const int type, const int k, const doublereal* const p); diff --git a/Cantera/src/transport/LiquidTransport.cpp b/Cantera/src/transport/LiquidTransport.cpp index 1b13824d0..903b1298f 100644 --- a/Cantera/src/transport/LiquidTransport.cpp +++ b/Cantera/src/transport/LiquidTransport.cpp @@ -39,10 +39,6 @@ namespace Cantera { m_iStateMF(-1), m_temp(-1.0), m_logt(0.0), - m_sqrt_t(-1.0), - m_t14(-1.0), - m_t32(-1.0), - m_sqrt_kbt(-1.0), m_press(-1.0), m_lambda(-1.0), m_viscmix(-1.0), @@ -67,10 +63,6 @@ namespace Cantera { m_iStateMF(-1), m_temp(-1.0), m_logt(0.0), - m_sqrt_t(-1.0), - m_t14(-1.0), - m_t32(-1.0), - m_sqrt_kbt(-1.0), m_press(-1.0), m_lambda(-1.0), m_viscmix(-1.0), @@ -100,36 +92,35 @@ namespace Cantera { m_tmin = right.m_tmin; m_tmax = right.m_tmax; m_mw = right.m_mw; - viscCoeffsVector_ = right.viscCoeffsVector_; - m_condcoeffs = right.m_condcoeffs; + m_visc_A = right.m_visc_A; + m_visc_logA = right.m_visc_logA; + m_visc_n = right.m_visc_n; + m_visc_Tact = right.m_visc_Tact; + m_visc_Eij = right.m_visc_Eij; + m_visc_Sij = right.m_visc_Sij; + m_thermCond_A = right.m_thermCond_A; + m_thermCond_n = right.m_thermCond_n; + m_thermCond_Tact = right.m_thermCond_Tact; + m_hydrodynamic_radius = right.m_hydrodynamic_radius; m_diffcoeffs = right.m_diffcoeffs; m_Grad_X = right.m_Grad_X; m_Grad_T = right.m_Grad_T; m_Grad_V = right.m_Grad_V; m_ck_Grad_mu = right.m_ck_Grad_mu; m_bdiff = right.m_bdiff; - viscSpecies_ = right.viscSpecies_; - m_sqvisc = right.m_sqvisc; - m_cond = right.m_cond; - m_polytempvec = right.m_polytempvec; + m_viscSpecies = right.m_viscSpecies; + m_logViscSpecies = right.m_logViscSpecies; + m_condSpecies = right.m_condSpecies; m_iStateMF = -1; m_molefracs = right.m_molefracs; m_concentrations = right.m_concentrations; m_chargeSpecies = right.m_chargeSpecies; m_DiffCoeff_StefMax = right.m_DiffCoeff_StefMax; viscosityModel_ = right.viscosityModel_; - m_phi = right.m_phi; - m_wratjk = right.m_wratjk; - m_wratkj1 = right.m_wratkj1; m_B = right.m_B; m_A = right.m_A; - m_eps = right.m_eps; m_temp = right.m_temp; m_logt = right.m_logt; - m_sqrt_t = right.m_sqrt_t; - m_t14 = right.m_t14; - m_t32 = right.m_t32; - m_sqrt_kbt = right.m_sqrt_kbt; m_press = right.m_press; m_flux = right.m_flux; m_lambda = right.m_lambda; @@ -173,32 +164,34 @@ namespace Cantera { copy(m_thermo->molecularWeights().begin(), m_thermo->molecularWeights().end(), m_mw.begin()); - // copy polynomials and parameters into local storage - viscCoeffsVector_ = tr.viscCoeffsVector_; - m_condcoeffs = tr.condcoeffs; + // copy parameters into local storage + m_visc_A = tr.visc_A ; + m_visc_n = tr.visc_n ; + m_visc_Tact = tr.visc_Tact ; + + //The following two are not yet filled in LiquidTransportParams + m_visc_Eij = tr.visc_Eij ; + m_visc_Sij = tr.visc_Sij ; + + //save logarithm of pre-exponential for easier computation + m_visc_logA.resize(m_nsp); + for ( int i = 0; i < m_nsp; i++ ) + m_visc_logA[i] = log( m_visc_A[i] ); + + m_thermCond_A = tr.thermCond_A ; + m_thermCond_n = tr.thermCond_n ; + m_thermCond_Tact = tr.thermCond_Tact ; + + m_hydrodynamic_radius = tr.hydroRadius ; + + //m_diffcoeffs = tr.diffcoeffs; - m_mode = tr.mode; - m_diam = tr.diam; - m_eps = tr.eps; + m_mode = tr.mode_; - m_phi.resize(m_nsp, m_nsp, 0.0); - - - m_wratjk.resize(m_nsp, m_nsp, 0.0); - m_wratkj1.resize(m_nsp, m_nsp, 0.0); - int j, k; - for (j = 0; j < m_nsp; j++) - for (k = j; k < m_nsp; k++) { - m_wratjk(j,k) = sqrt(m_mw[j]/m_mw[k]); - m_wratjk(k,j) = sqrt(m_wratjk(j,k)); - m_wratkj1(j,k) = sqrt(1.0 + m_mw[k]/m_mw[j]); - } - - m_polytempvec.resize(5); - viscSpecies_.resize(m_nsp); - m_sqvisc.resize(m_nsp); - m_cond.resize(m_nsp); + m_viscSpecies.resize(m_nsp); + m_logViscSpecies.resize(m_nsp); + m_condSpecies.resize(m_nsp); m_bdiff.resize(m_nsp, m_nsp); m_molefracs.resize(m_nsp); @@ -249,7 +242,7 @@ namespace Cantera { if (m_visc_mix_ok) return m_viscmix; - // update viscSpecies_[] and m_phi[] if necessary + // update m_viscSpecies[] if necessary if (!m_visc_temp_ok) { updateViscosity_temp(); } @@ -258,16 +251,27 @@ namespace Cantera { updateViscosities_conc(); } + /* We still need to implement interaction parameters */ + /* This constant viscosity model has no input */ if (viscosityModel_ == LVISC_CONSTANT) { - return m_viscmix; - } else if (viscosityModel_ == LVISC_MIXTUREAVG) { - m_viscmix = dot_product(viscSpecies_, m_molefracs); - } else if (viscosityModel_ == LVISC_WILKES) { - multiply(m_phi, DATA_PTR(m_molefracs), DATA_PTR(m_spwork)); - m_viscmix = 0.0; - for (int k = 0; k < m_nsp; k++) { - m_viscmix += m_molefracs[k] * viscSpecies_[k]/m_spwork[k]; - } + + err("constant viscosity not implemented for LiquidTransport."); + //return m_viscmix; + + } else if (viscosityModel_ == LVISC_AVG_ENERGIES) { + + m_viscmix = exp( dot_product(m_logViscSpecies, m_molefracs) ); + + } else if (viscosityModel_ == LVISC_INTERACTION) { + + // log_visc_mix = sum_i (X_i log_visc_i) + sum_i sum_j X_i X_j G_ij + double interaction = dot_product(m_logViscSpecies, m_molefracs); + for ( int i = 0; i < m_nsp; i++ ) + for ( int j = 0; j < i; j++ ) + interaction += m_molefracs[i] * m_molefracs[j] + * ( m_visc_Sij(i,j) + m_visc_Eij(i,j) / m_temp ); + m_viscmix = exp( interaction ); + } return m_viscmix; @@ -278,7 +282,7 @@ namespace Cantera { if (!m_visc_temp_ok) { updateViscosity_temp(); } - copy(viscSpecies_.begin(), viscSpecies_.end(), visc); + copy(m_viscSpecies.begin(), m_viscSpecies.end(), visc); } @@ -301,41 +305,79 @@ namespace Cantera { d[ld*j + i] = rp * m_bdiff(i,j); } } - - + //================================================================================================ + // Get the electrical Mobilities (m^2/V/s). + /* + * This function returns the mobilities. In some formulations + * this is equal to the normal mobility multiplied by faraday's constant. + * + * Frequently, but not always, the mobility is calculated from the + * diffusion coefficient using the Einstein relation + * + * \f[ + * \mu^e_k = \frac{F D_k}{R T} + * \f] + * + * @param mobil_e Returns the mobilities of + * the species in array \c mobil_e. The array must be + * dimensioned at least as large as the number of species. + */ void LiquidTransport::getMobilities(doublereal* const mobil) { - // this needs to be checked out. int k; getMixDiffCoeffs(DATA_PTR(m_spwork)); doublereal c1 = ElectronCharge / (Boltzmann * m_temp); for (k = 0; k < m_nsp; k++) { - mobil[k] = c1 * m_spwork[k] * m_thermo->charge(k); + mobil[k] = c1 * m_spwork[k]; } } - - + //================================================================================================ + //! Get the fluid mobilities (s kmol/kg). + /*! + * This function returns the fluid mobilities. Usually, you have + * to multiply Faraday's constant into the resulting expression + * to general a species flux expression. + * + * Frequently, but not always, the mobility is calculated from the + * diffusion coefficient using the Einstein relation + * + * \f[ + * \mu^f_k = \frac{D_k}{R T} + * \f] + * + * + * @param mobil_f Returns the mobilities of + * the species in array \c mobil. The array must be + * dimensioned at least as large as the number of species. + */ + void LiquidTransport::getFluidMobilities(doublereal* const mobil_f) { + getMixDiffCoeffs(DATA_PTR(m_spwork)); + doublereal c1 = 1.0 / (GasConstant * m_temp); + for (int k = 0; k < m_nsp; k++) { + mobil_f[k] = c1 * m_spwork[k]; + } + } + //================================================================================================ void LiquidTransport::set_Grad_V(const doublereal* const grad_V) { for (int a = 0; a < m_nDim; a++) { m_Grad_V[a] = grad_V[a]; } } - + //================================================================================================ void LiquidTransport::set_Grad_T(const doublereal* const grad_T) { for (int a = 0; a < m_nDim; a++) { m_Grad_T[a] = grad_T[a]; } } - - void LiquidTransport::set_Grad_X(const doublereal* const grad_X) { - int itop = m_nDim * m_nsp; - for (int i = 0; i < itop; i++) { - m_Grad_X[i] = grad_X[i]; - } - update_Grad_lnAC(); - } - - + //================================================================================================ + void LiquidTransport::set_Grad_X(const doublereal* const grad_X) { + int itop = m_nDim * m_nsp; + for (int i = 0; i < itop; i++) { + m_Grad_X[i] = grad_X[i]; + } + update_Grad_lnAC(); + } + //================================================================================================ /****************** thermal conductivity **********************/ /* @@ -356,8 +398,8 @@ namespace Cantera { if (!m_cond_mix_ok) { doublereal sum1 = 0.0, sum2 = 0.0; for (int k = 0; k < m_nsp; k++) { - sum1 += m_molefracs[k] * m_cond[k]; - sum2 += m_molefracs[k] / m_cond[k]; + sum1 += m_molefracs[k] * m_condSpecies[k]; + sum2 += m_molefracs[k] / m_condSpecies[k]; } m_lambda = 0.5*(sum1 + 1.0/sum2); m_cond_mix_ok = true; @@ -439,35 +481,6 @@ namespace Cantera { } } - void LiquidTransport::getSpeciesDiffusiveMassFluxes(doublereal* const fluxes) { - int n, k; - - update_temp(); - update_conc(); - - - getMixDiffCoeffs(DATA_PTR(m_spwork)); - - - const array_fp& mw = m_thermo->molecularWeights(); - const doublereal* const y = m_thermo->massFractions(); - const doublereal rhon = m_thermo->molarDensity(); - // Unroll wrt ndim - vector_fp sum(m_nDim,0.0); - for (n = 0; n < m_nDim; n++) { - for (k = 0; k < m_nsp; k++) { - fluxes[n*m_nsp + k] = -rhon * mw[k] * m_spwork[k] * m_Grad_X[n*m_nsp + k]; - sum[n] += fluxes[n*m_nsp + k]; - } - } - // add correction flux to enforce sum to zero - for (n = 0; n < m_nDim; n++) { - for (k = 0; k < m_nsp; k++) { - fluxes[n*m_nsp + k] -= y[k]*sum[n]; - } - } - } - /** * Mixture-averaged diffusion coefficients [m^2/s]. * @@ -542,18 +555,6 @@ namespace Cantera { m_temp = t; m_logt = log(m_temp); m_kbt = Boltzmann * m_temp; - m_sqrt_t = sqrt(m_temp); - m_t14 = sqrt(m_sqrt_t); - m_t32 = m_temp * m_sqrt_t; - m_sqrt_kbt = sqrt(Boltzmann*m_temp); - - // compute powers of log(T) - // -> may move this - m_polytempvec[0] = 1.0; - m_polytempvec[1] = m_logt; - m_polytempvec[2] = m_logt*m_logt; - m_polytempvec[3] = m_logt*m_logt*m_logt; - m_polytempvec[4] = m_logt*m_logt*m_logt*m_logt; // temperature has changed so temp flags are flipped m_visc_temp_ok = false; @@ -689,18 +690,20 @@ namespace Cantera { */ void LiquidTransport::updateCond_temp() { - int k; + + /* if (m_mode == CK_Mode) { for (k = 0; k < m_nsp; k++) { - m_cond[k] = exp(dot4(m_polytempvec, m_condcoeffs[k])); + m_condSpecies[k] = exp(m_condcoeffs[k]); } } else { for (k = 0; k < m_nsp; k++) { - m_cond[k] = m_sqrt_t * dot5(m_polytempvec, m_condcoeffs[k]); + m_condSpecies[k] = m_sqrt_t * m_condcoeffs[k]; } } m_cond_temp_ok = true; m_cond_mix_ok = false; + */ } @@ -711,12 +714,12 @@ namespace Cantera { void LiquidTransport::updateDiff_temp() { // evaluate binary diffusion coefficients at unit pressure - int i,j; - int ic = 0; + + /* if (m_mode == CK_Mode) { for (i = 0; i < m_nsp; i++) { for (j = i; j < m_nsp; j++) { - m_bdiff(i,j) = exp(dot4(m_polytempvec, m_diffcoeffs[ic])); + m_bdiff(i,j) = exp(m_diffcoeffs[ic]); m_bdiff(j,i) = m_bdiff(i,j); ic++; } @@ -725,8 +728,7 @@ namespace Cantera { else { for (i = 0; i < m_nsp; i++) { for (j = i; j < m_nsp; j++) { - m_bdiff(i,j) = m_temp * m_sqrt_t*dot5(m_polytempvec, - m_diffcoeffs[ic]); + m_bdiff(i,j) = m_temp * m_sqrt_t*m_diffcoeffs[ic]; m_bdiff(j,i) = m_bdiff(i,j); ic++; } @@ -735,6 +737,7 @@ namespace Cantera { m_diff_temp_ok = true; m_diff_mix_ok = false; + */ } @@ -754,38 +757,16 @@ namespace Cantera { */ void LiquidTransport::updateViscosity_temp() { int k; - doublereal vratiokj, wratiojk, factor1; - if (m_mode == CK_Mode) { - for (k = 0; k < m_nsp; k++) { - viscSpecies_[k] = exp(dot4(m_polytempvec, viscCoeffsVector_[k])); - m_sqvisc[k] = sqrt(viscSpecies_[k]); - } + for (k = 0; k < m_nsp; k++) { + m_logViscSpecies[k] = m_visc_logA[k] + m_visc_n[k] * m_logt + + m_visc_Tact[k] / m_temp ; + m_viscSpecies[k] = exp( m_logViscSpecies[k] ); } - else { - for (k = 0; k < m_nsp; k++) { - // the polynomial fit is done for sqrt(visc/sqrt(T)) - m_sqvisc[k] = m_t14*dot5(m_polytempvec, viscCoeffsVector_[k]); - viscSpecies_[k] = (m_sqvisc[k]*m_sqvisc[k]); - } - } - - // see Eq. (9-5.15) of Reid, Prausnitz, and Poling - int j; - for (j = 0; j < m_nsp; j++) { - for (k = j; k < m_nsp; k++) { - vratiokj = viscSpecies_[k]/viscSpecies_[j]; - wratiojk = m_mw[j]/m_mw[k]; - - // Note that m_wratjk(k,j) holds the square root of - // m_wratjk(j,k)! - factor1 = 1.0 + (m_sqvisc[k]/m_sqvisc[j]) * m_wratjk(k,j); - m_phi(k,j) = factor1*factor1 / - (SqrtEight * m_wratkj1(j,k)); - m_phi(j,k) = m_phi(k,j)/(vratiokj * wratiojk); - } - } - + //for (k = 0; k < m_nsp; k++) { + //m_viscSpecies[k] = m_visc_A[k] * exp( m_visc_n[k] * m_logt + // + m_visc_Tact[k] / m_temp ); + //} m_visc_temp_ok = true; m_visc_mix_ok = false; } @@ -953,4 +934,20 @@ namespace Cantera { } } } + + + /** + * Throw an exception if this method is invoked. + * This probably indicates something is not yet implemented. + */ + doublereal LiquidTransport::err(std::string msg) const { + throw CanteraError("Liquid Transport Class", + "\n\n\n**** Method "+ msg +" not implemented in model " + + int2str(model()) + " ****\n" + "(Did you forget to specify a transport model?)\n\n\n"); + + return 0.0; + } + + } diff --git a/Cantera/src/transport/LiquidTransport.h b/Cantera/src/transport/LiquidTransport.h index 1b131d764..c3df40c3f 100644 --- a/Cantera/src/transport/LiquidTransport.h +++ b/Cantera/src/transport/LiquidTransport.h @@ -25,13 +25,14 @@ using namespace std; // Cantera includes #include "TransportBase.h" #include "DenseMatrix.h" +#include "TransportParams.h" #include "LiquidTransportParams.h" namespace Cantera { const int LVISC_CONSTANT = 0; - const int LVISC_WILKES = 1; - const int LVISC_MIXTUREAVG = 2; + const int LVISC_INTERACTION = 1; + const int LVISC_AVG_ENERGIES = 2; const int LDIFF_MIXDIFF_UNCORRECTED = 0; const int LDIFF_MIXDIFF_FLUXCORRECTED = 1; @@ -39,7 +40,7 @@ namespace Cantera { - class TransportParams; + class LiquidTransportParams; //! Class LiquidTransport implements mixture-averaged transport @@ -135,12 +136,21 @@ namespace Cantera { class LiquidTransport : public Transport { public: - //! default constructor + //! Default constructor. + /*! + * This requires call to initLiquid(LiquidTransportParams& tr) + * after filling LiquidTransportParams to complete instantiation. + * The filling of LiquidTransportParams is currently carried out + * in the TransportFactory class, but might be moved at some point. + * + * @param thermo ThermoPhase object holding species information. + * @param ndim Number of spatial dimensions. + */ LiquidTransport(thermo_t* thermo = 0, int ndim = 1); //!Copy Constructor for the %LiquidThermo object. /*! - * @param right ThermoPhase to be copied + * @param right %LiquidTransport to be copied */ LiquidTransport(const LiquidTransport &right); @@ -148,8 +158,8 @@ namespace Cantera { /*! * This is NOT a virtual function. * - * @param right Reference to %ThermoPhase object to be copied into the - * current one. + * @param right Reference to %LiquidTransport object to be copied + * into the current one. */ LiquidTransport& operator=(const LiquidTransport& right); @@ -169,8 +179,21 @@ namespace Cantera { //! virtual destructor virtual ~LiquidTransport() {} + //! Initialize the transport object + /*! + * Here we change all of the internal dimensions to be sufficient. + * We get the object ready to do property evaluations. + * + * @param tr Transport parameters for all of the species + * in the phase. + */ + virtual bool initLiquid(LiquidTransportParams& tr); + + friend class TransportFactory; + + //! Return the model id for this transport parameterization - virtual int model() { + virtual int model() const { return cLiquidTransport; } @@ -178,17 +201,14 @@ namespace Cantera { //! Returns the viscosity of the solution /*! - * The viscosity is computed using the Wilke mixture rule. + * The viscosity is computed using mixture averaging plus + * any information on interaction parameters * \f[ - * \mu = \sum_k \frac{\mu_k X_k}{\sum_j \Phi_{k,j} X_j}. + * \mu = \sum_k {\mu_k X_k} {\sum_j \sum_k {G_{j,k} X_k X_j} }. * \f] * Here \f$ \mu_k \f$ is the viscosity of pure species \e k, - * 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{8}\sqrt{1 + M_k/M_j}} - * \f] + * and \f$ G_{k,j} \f$ is the interaction parameter. + * @see updateViscosity_T(); * * Controlling update boolean m_viscmix_ok @@ -197,25 +217,11 @@ namespace Cantera { //! Returns the pure species viscosities /*! - * - * + * The pure species viscosities are to be given in an Arrhenius + * form in accordance with activated-jump-process dominated transport. */ virtual void getSpeciesViscosities(doublereal* const visc); - virtual void getThermalDiffCoeffs(doublereal* const dt); - - //! Return the thermal conductivity of the solution - /*! - * The thermal conductivity is computed from the following mixture rule: - * \f[ - * \lambda = 0.5 \left( \sum_k X_k \lambda_k - * + \frac{1}{\sum_k X_k/\lambda_k}\right) - * \f] - * - * Controlling update boolean = m_condmix_ok - */ - virtual doublereal thermalConductivity(); - //! Returns the binary diffusion coefficients /*! * @param ld @@ -231,11 +237,56 @@ namespace Cantera { virtual void getMixDiffCoeffs(doublereal* const d); - //! Get the Mobilities + virtual void getThermalDiffCoeffs(doublereal* const dt); + + //! Return the thermal conductivity of the solution /*! - * @param mobil + * The thermal conductivity is computed from the following mixture rule: + * \f[ + * \lambda = 0.5 \left( \sum_k X_k \lambda_k + * + \frac{1}{\sum_k X_k/\lambda_k}\right) + * \f] + * + * Controlling update boolean = m_condmix_ok */ - virtual void getMobilities(doublereal* const mobil); + virtual doublereal thermalConductivity(); + + //! Get the Electrical mobilities (m^2/V/s). + /*! + * This function returns the mobilities. In some formulations + * this is equal to the normal mobility multiplied by faraday's constant. + * + * The mobility is calculated from the + * diffusion coefficient using the Einstein relation + * + * \f[ + * \mu^e_k = \frac{F D_k}{R T} + * \f] + * + * @param mobil_e Returns the electrical mobilities of + * the species in array \c mobil_e. The array must be + * dimensioned at least as large as the number of species. + */ + virtual void getMobilities(doublereal* const mobil_e); + + //! Get the fluid mobilities (s kmol/kg). + /*! + * This function returns the fluid mobilities. Usually, you have + * to multiply Faraday's constant into the resulting expression + * to general a species flux expression. + * + * The mobility is calculated from the + * diffusion coefficient using the Einstein relation + * + * \f[ + * \mu^f_k = \frac{D_k}{R T} + * \f] + * + * @param mobil_f Returns the fluid mobilities of + * the species in array \c mobil_f. The array must be + * dimensioned at least as large as the number of species. + */ + virtual void getFluidMobilities(doublereal* const mobil_f); //! Specify the value of the gradient of the voltage /*! @@ -306,19 +357,6 @@ namespace Cantera { int ldx, const doublereal* grad_X, int ldf, doublereal* fluxes); - //! Return the species diffusive mass fluxes - /*! - * - * - * - * @param ndim The number of spatial dimensions (1, 2, or 3). - * @param grad_T The temperature gradient (ignored in this model). - * @param ldx Leading dimension of the grad_X array. - * The diffusive mass flux of species \e k is computed from - * - * - */ - virtual void getSpeciesDiffusiveMassFluxes(doublereal* const fluxes); /** * @param ndim The number of spatial dimensions (1, 2, or 3). @@ -331,20 +369,6 @@ namespace Cantera { virtual void getSpeciesFluxesExt(int ldf, doublereal* fluxes); - //! Initialize the transport object - /*! - * Here we change all of the internal dimensions to be sufficient. - * We get the object ready to do property evaluations. - * - * @param tr Transport parameters for all of the species - * in the phase. - */ - virtual bool initLiquid(LiquidTransportParams& tr); - - friend class TransportFactory; - - - //! Solve the stefan_maxell equations for the diffusive fluxes. void stefan_maxwell_solve(); @@ -368,19 +392,34 @@ namespace Cantera { */ vector_fp m_mw; - //! Polynomial coefficients of the viscosity - /*! - * These express the temperature dependendence of the pures - * species viscosities. - */ - std::vector viscCoeffsVector_; + //! Pure species viscosities in Arrhenius temperature-dependent form. + vector_fp m_visc_A; + vector_fp m_visc_logA; //logarithm of coefficient + vector_fp m_visc_n; + vector_fp m_visc_Tact; - //! Polynomial coefficients of the conductivities - /*! - * These express the temperature dependendence of the pures - * species conductivities + //! Molecular interaction energies associated with viscosity + /** + * These multiply the viscosity according to + * \f[ exp( \sum_{i} \sum_{j} X_i X_j E_{i,j} / T \f]. */ - vector m_condcoeffs; + DenseMatrix m_visc_Eij; + + //! Molecular interaction entropies associated with viscosity + /** + * These multiply the viscosity according to + * \f[ exp( \sum_{i} \sum{j} X_i X_j S_{i,j} \f]. + */ + DenseMatrix m_visc_Sij; + + //! Pure species thermal conductivities in Arrhenius temperature-dependent form. + vector_fp m_thermCond_A; + vector_fp m_thermCond_n; + vector_fp m_thermCond_Tact; + + //! Species hydrodynamic radius + vector_fp m_hydrodynamic_radius; + //! Polynomial coefficients of the binary diffusion coefficients /*! @@ -389,7 +428,6 @@ namespace Cantera { * added. */ vector m_diffcoeffs; - //! Internal value of the gradient of the mole fraction vector /*! @@ -472,9 +510,9 @@ namespace Cantera { */ DenseMatrix m_bdiff; - //! Species viscosities + //! Species viscosities and their logarithm /*! - * Viscosity of the species + * Viscosity of the species and its logarithm * Length = number of species * * Depends on the temperature. We have set the pressure dependence @@ -482,19 +520,8 @@ namespace Cantera { * * controlling update boolean -> m_visc_temp_ok */ - vector_fp viscSpecies_; - - //! Sqrt of the species viscosities - /*! - * The sqrt(visc) is used in the mixing formulas - * Length = m_nsp - * - * Depends on the temperature and perhaps pressure, but - * not the species concentrations - * - * controlling update boolean m_visc_temp_ok - */ - vector_fp m_sqvisc; + vector_fp m_viscSpecies; + vector_fp m_logViscSpecies; //! Internal value of the species individual thermal conductivities /*! @@ -505,10 +532,7 @@ namespace Cantera { * * controlling update boolean -> m_cond_temp_ok */ - vector_fp m_cond; - - //! Polynomials of the log of the temperature - vector_fp m_polytempvec; + vector_fp m_condSpecies; //! State of the mole fraction vector. int m_iStateMF; @@ -592,31 +616,12 @@ namespace Cantera { */ int viscosityModel_; - //! viscosity weighting functions - DenseMatrix m_phi; - - //! Matrix of the ratios of the species molecular weights - /*! - * m_wratjk(i,j) = (m_mw[j]/m_mw[k])**0.25 - */ - DenseMatrix m_wratjk; - - //! Matrix of the ratios of the species molecular weights - /*! - * m_wratkj1(i,j) = (1.0 + m_mw[k]/m_mw[j])**0.5 - */ - DenseMatrix m_wratkj1; - //! RHS to the stefan-maxwell equation DenseMatrix m_B; //! Matrix for the stefan maxwell equation. DenseMatrix m_A; - //! Internal storage for the species LJ well depth - vector_fp m_eps; - - //! Current Temperature -> locally storred /*! * This is used to test whether new temperature computations @@ -630,21 +635,6 @@ namespace Cantera { //! Current value of kT doublereal m_kbt; - //! Current Temperature **0.5 - doublereal m_sqrt_t; - - //! Current Temperature **0.25 - doublereal m_t14; - - //! Current Temperature **1.5 - doublereal m_t32; - - //! Current temperature function - /*! - * This is equal to sqrt(Boltzmann * T) - */ - doublereal m_sqrt_kbt; - //! Current value of the pressure doublereal m_press; @@ -745,6 +735,17 @@ namespace Cantera { * Either 1, 2, or 3 */ int m_nDim; + + private: + + //! Throw an exception if this method is invoked. + /*! + * This probably indicates something is not yet implemented. + * + * @pram msg Indicates the member function which is not implemented + */ + doublereal err(std::string msg) const; + }; } #endif diff --git a/Cantera/src/transport/LiquidTransportData.h b/Cantera/src/transport/LiquidTransportData.h new file mode 100644 index 000000000..a664bedc1 --- /dev/null +++ b/Cantera/src/transport/LiquidTransportData.h @@ -0,0 +1,80 @@ +/** + * @file TransportFactory.h + * Header file defining class TransportFactory + * (see \link Cantera::TransportFactory TransportFactory\endlink) + */ +/* + * $Author: hkmoffa $ + * $Date: 2008/12/24 18:19:01 $ + * $Revision: 1.14 $ + * + * + * + */ + +#ifndef CT_LIQUIDTRANSPORTDATA_H +#define CT_LIQUIDTRANSPORTDATA_H + + +// STL includes +#include +#include +#include +#include + + + +// Cantera includes +#include "ct_defs.h" +#include "TransportBase.h" +#include "FactoryBase.h" + + +namespace Cantera { + + enum LiquidTR_Model { + LTR_MODEL_NOTSET=-1, + LTR_MODEL_CONSTANT, + LTR_MODEL_ARRHENIUS, + LTR_MODEL_COEFF + }; + + class LiquidTransportData { + + public: + + LiquidTransportData() : + speciesName("-"), + model_hydroradius(LTR_MODEL_NOTSET), + hydroradius(-1.0), + model_viscosity(LTR_MODEL_NOTSET), + model_thermalCond(LTR_MODEL_NOTSET), + model_speciesDiffusivity(LTR_MODEL_NOTSET) + { + } + + std::string speciesName; + + //! Model type for the hydroradius + LiquidTR_Model model_hydroradius; + + //! Actual value of the hydroradius + doublereal hydroradius; + + //! Model type for the hydroradius + LiquidTR_Model model_viscosity; + vector_fp viscCoeffs; + + //! Model type for the hydroradius + LiquidTR_Model model_thermalCond; + + vector_fp thermalCondCoeffs; + + //! Model type for the hydroradius + LiquidTR_Model model_speciesDiffusivity; + + vector_fp speciesDiffusivityCoeffs; + }; + +} +#endif diff --git a/Cantera/src/transport/LiquidTransportParams.h b/Cantera/src/transport/LiquidTransportParams.h index e55f0bf58..9b0d98c94 100644 --- a/Cantera/src/transport/LiquidTransportParams.h +++ b/Cantera/src/transport/LiquidTransportParams.h @@ -5,27 +5,62 @@ #include "ct_defs.h" #include "TransportBase.h" +#include "TransportParams.h" +#include "LiquidTransportData.h" #include "xml.h" #include "XML_Writer.h" namespace Cantera { /** - * - * Holds transport data. Used by TransportFactory. - * + * Holds transport model parameters relevant to transport in + * liquids for which activated jump processes limit transport + * (giving Arrhenius type transport properties). + * Used by TransportFactory. */ - class LiquidTransportParams { + class LiquidTransportParams :public TransportParams { public: - LiquidTransportParams() : thermo(0), xml(0) {} - virtual ~LiquidTransportParams(); - int nsp; + LiquidTransportParams() {} + ~LiquidTransportParams() {} - // phase_t* mix; - thermo_t* thermo; - vector_fp mw; + + //section for liquid transport properties + + //Arrhenius parameters for transport coefficients: + + //!Arrhenius pre-exponential parameter for viscosity. + vector_fp visc_A; + //!Temperature exponent for viscosity. + vector_fp visc_n; + //!Arrhenius activation temperature for viscosity. + vector_fp visc_Tact; + + //!Arrhenius pre-exponential parameter for thermal conductivity. + vector_fp thermCond_A; + //!Temperature exponent for thermal conductivity. + vector_fp thermCond_n; + //!Arrhenius activation temperature for thermal conductivity. + vector_fp thermCond_Tact; + + //! Energies of molecular interaction associated with viscosity. + /** + * These multiply the mixture viscosity by + * \f[ \exp( \sum_{i} \sum_{j} X_i X_j ( S_{i,j} + E_{i,j} / T ) ) \f]. + * + * The overall formula for the logarithm of the mixture viscosity is + * + * \f[ \ln \eta_{mix} = \sum_i X_i \ln \eta_i + * + \sum_i \sum_j X_i X_j ( S_{i,j} + E_{i,j} / T ) \f]. + */ + DenseMatrix visc_Eij; + + //! Entropies of molecular interaction associated with viscosity. + DenseMatrix visc_Sij; + + //Hydrodynamic radius of transported molecule + vector_fp hydroRadius; //! Coefficients for the limiting conductivity of ions //! in solution: A_k @@ -47,26 +82,8 @@ namespace Cantera { vector_fp B_k_cond; - // polynomial fits - std::vector viscCoeffsVector_; - std::vector condcoeffs; - std::vector diffcoeffs ; + std::vector LTData; - - std::vector polar; - //vector_fp alpha; - vector_fp fitlist; - vector_fp eps; - vector_fp sigma; - DenseMatrix reducedMass; - DenseMatrix diam; - DenseMatrix epsilon; - DenseMatrix dipole; - DenseMatrix delta; - doublereal tmax, tmin; - int mode; - XML_Writer* xml; - int log_level; }; } diff --git a/Cantera/src/transport/Makefile.in b/Cantera/src/transport/Makefile.in index ad6d76b15..9f52f08e1 100644 --- a/Cantera/src/transport/Makefile.in +++ b/Cantera/src/transport/Makefile.in @@ -35,11 +35,13 @@ CXX_FLAGS = @CXXFLAGS@ $(CXX_OPT) $(PIC_FLAG) $(DEBUG_FLAG) # Base Transport Object Files TRAN_OBJ = TransportFactory.o MultiTransport.o MixTransport.o MMCollisionInt.o \ - SolidTransport.o DustyGasTransport.o TransportBase.o WaterTransport.o + SolidTransport.o DustyGasTransport.o TransportBase.o WaterTransport.o \ + SimpleTransport.o TRAN_H = TransportFactory.h MultiTransport.h MixTransport.h \ MMCollisionInt.h SolidTransport.h DustyGasTransport.h \ - TransportBase.h L_matrix.h TransportParams.h WaterTransport.h + TransportBase.h L_matrix.h TransportParams.h WaterTransport.h \ + SimpleTransport.h LiquidTransportData.h ifeq ($(do_electro),1) do_issp = 1 diff --git a/Cantera/src/transport/MixTransport.cpp b/Cantera/src/transport/MixTransport.cpp index c7fe8ef8a..2dd2da5da 100755 --- a/Cantera/src/transport/MixTransport.cpp +++ b/Cantera/src/transport/MixTransport.cpp @@ -50,7 +50,7 @@ namespace Cantera { } - bool MixTransport::init(TransportParams& tr) { + bool MixTransport::initGas( GasTransportParams& tr ) { // constant substance attributes m_thermo = tr.thermo; @@ -72,7 +72,7 @@ namespace Cantera { m_zrot = tr.zrot; m_crot = tr.crot; m_epsilon = tr.epsilon; - m_mode = tr.mode; + m_mode = tr.mode_; m_diam = tr.diam; m_eps = tr.eps; m_alpha = tr.alpha; diff --git a/Cantera/src/transport/MixTransport.h b/Cantera/src/transport/MixTransport.h old mode 100755 new mode 100644 index 5fc8d4b96..18a781577 --- a/Cantera/src/transport/MixTransport.h +++ b/Cantera/src/transport/MixTransport.h @@ -37,7 +37,7 @@ using namespace std; namespace Cantera { - class TransportParams; + class GasTransportParams; /** * Class MixTransport implements mixture-averaged transport @@ -51,7 +51,7 @@ namespace Cantera { virtual ~MixTransport() {} - virtual int model() { return cMixtureAveraged; } + virtual int model() const { return cMixtureAveraged; } //! Viscosity of the mixture /*! @@ -120,7 +120,7 @@ namespace Cantera { * @param tr Transport parameters for all of the species * in the phase. */ - virtual bool init(TransportParams& tr); + virtual bool initGas( GasTransportParams& tr ); friend class TransportFactory; diff --git a/Cantera/src/transport/MultiTransport.cpp b/Cantera/src/transport/MultiTransport.cpp index 20bf2022e..1f9e443c0 100755 --- a/Cantera/src/transport/MultiTransport.cpp +++ b/Cantera/src/transport/MultiTransport.cpp @@ -142,7 +142,7 @@ namespace Cantera { } - bool MultiTransport::init(TransportParams& tr) { + bool MultiTransport::initGas( GasTransportParams& tr ) { // constant mixture attributes //m_phase = tr.mix; @@ -167,7 +167,7 @@ namespace Cantera { m_zrot = tr.zrot; m_crot = tr.crot; m_epsilon = tr.epsilon; - m_mode = tr.mode; + m_mode = tr.mode_; m_diam = tr.diam; m_eps = tr.eps; m_alpha = tr.alpha; diff --git a/Cantera/src/transport/MultiTransport.h b/Cantera/src/transport/MultiTransport.h old mode 100755 new mode 100644 index 78a0f998a..5979b8ce9 --- a/Cantera/src/transport/MultiTransport.h +++ b/Cantera/src/transport/MultiTransport.h @@ -40,7 +40,7 @@ namespace Cantera { TRANSOLVE_LU }; - class TransportParams; + class GasTransportParams; ///////////////////////////////////////////////////////////// @@ -85,7 +85,7 @@ namespace Cantera { virtual ~MultiTransport(); // overloaded base class methods - virtual int model() { + virtual int model() const { if (m_mode == CK_Mode) return CK_Multicomponent; else @@ -162,7 +162,7 @@ namespace Cantera { /** * @internal */ - virtual bool init(TransportParams& tr); + virtual bool initGas( GasTransportParams& tr ); /** diff --git a/Cantera/src/transport/SimpleTransport.cpp b/Cantera/src/transport/SimpleTransport.cpp new file mode 100644 index 000000000..eafaa02fe --- /dev/null +++ b/Cantera/src/transport/SimpleTransport.cpp @@ -0,0 +1,758 @@ +/** + * @file SimpleTransport.cpp + * Simple mostly constant transport properties + */ +/* + * $Revision: 1.10 $ + * $Date: 2009/03/24 20:44:30 $ + */ + +#include "ThermoPhase.h" +#include "SimpleTransport.h" + +#include "utilities.h" +#include "LiquidTransportParams.h" +#include "TransportFactory.h" + +#include "ctlapack.h" + +#include +using namespace std; + +/** + * Mole fractions below MIN_X will be set to MIN_X when computing + * transport properties. + */ +#define MIN_X 1.e-14 + + +namespace Cantera { + //================================================================================================ + SimpleTransport::SimpleTransport(thermo_t* thermo, int ndim) : + Transport(thermo, ndim), + m_nsp(0), + tempDepType_(0), + compositionDepType_(0), + useHydroRadius_(false), + doMigration_(0), + m_tmin(-1.0), + m_tmax(100000.), + m_iStateMF(-1), + concTot_(0.0), + m_temp(-1.0), + m_press(-1.0), + m_lambda(-1.0), + m_viscmix(-1.0), + m_visc_mix_ok(false), + m_visc_temp_ok(false), + m_diff_mix_ok(false), + m_diff_temp_ok(false), + m_cond_temp_ok(false), + m_cond_mix_ok(false), + m_nDim(1) + { + } + //================================================================================================ + SimpleTransport::SimpleTransport(const SimpleTransport &right) : + Transport(), + m_nsp(0), + tempDepType_(0), + compositionDepType_(0), + useHydroRadius_(false), + doMigration_(0), + m_tmin(-1.0), + m_tmax(100000.), + m_iStateMF(-1), + m_temp(-1.0), + m_press(-1.0), + m_lambda(-1.0), + m_viscmix(-1.0), + m_visc_mix_ok(false), + m_visc_temp_ok(false), + m_diff_mix_ok(false), + m_diff_temp_ok(false), + m_cond_temp_ok(false), + m_cond_mix_ok(false), + m_nDim(1) + { + /* + * Use the assignment operator to do the brunt + * of the work for the copy construtor. + */ + *this = right; + } + //================================================================================================ + SimpleTransport& SimpleTransport::operator=(const SimpleTransport& right) { + if (&right != this) { + return *this; + } + Transport::operator=(right); + + m_nsp = right.m_nsp; + tempDepType_ = right.tempDepType_; + compositionDepType_ = right.compositionDepType_; + useHydroRadius_ = right.useHydroRadius_; + doMigration_ = right.doMigration_; + m_tmin = right.m_tmin; + m_tmax = right.m_tmax; + m_mw = right.m_mw; + + m_coeffVisc_Ns = right.m_coeffVisc_Ns; + m_coeffLambda_Ns = right.m_coeffLambda_Ns; + m_coeffDiff_Ns = right.m_coeffDiff_Ns; + + m_Grad_X = right.m_Grad_X; + m_Grad_T = right.m_Grad_T; + m_Grad_P = right.m_Grad_P; + m_Grad_V = right.m_Grad_V; + + m_diffSpecies = right.m_diffSpecies; + m_viscSpecies = right.m_viscSpecies; + m_condSpecies = right.m_condSpecies; + m_iStateMF = -1; + m_molefracs = right.m_molefracs; + m_concentrations = right.m_concentrations; + concTot_ = right.concTot_; + meanMolecularWeight_ = right.meanMolecularWeight_; + dens_ = right.dens_; + m_chargeSpecies = right.m_chargeSpecies; + + m_temp = right.m_temp; + m_press = right.m_press; + m_lambda = right.m_lambda; + m_viscmix = right.m_viscmix; + m_spwork = right.m_spwork; + m_visc_mix_ok = false; + m_visc_temp_ok = false; + m_diff_mix_ok = false; + m_diff_temp_ok = false; + m_cond_temp_ok = false; + m_cond_mix_ok = false; + m_nDim = right.m_nDim; + + return *this; + } + + //================================================================================================ + Transport *SimpleTransport::duplMyselfAsTransport() const { + SimpleTransport* tr = new SimpleTransport(*this); + return (dynamic_cast(tr)); + } + //================================================================================================ + // Initialize the object + /* + * This is where we dimension everything. + */ + bool SimpleTransport::initLiquid(LiquidTransportParams& tr) { + int k; + // constant substance attributes + m_thermo = tr.thermo; + m_nsp = m_thermo->nSpecies(); + m_tmin = m_thermo->minTemp(); + m_tmax = m_thermo->maxTemp(); + + // make a local copy of the molecular weights + m_mw.resize(m_nsp); + copy(m_thermo->molecularWeights().begin(), + m_thermo->molecularWeights().end(), m_mw.begin()); + + /* + * Get the input Viscosities + */ + m_viscSpecies.resize(m_nsp); + m_coeffVisc_Ns.clear(); + m_coeffVisc_Ns.resize(m_nsp); + + Cantera::LiquidTransportData <d0 = tr.LTData[0]; + LiquidTR_Model vm0 = ltd0.model_viscosity; + if (vm0 == LTR_MODEL_CONSTANT) { + tempDepType_ = 0; + } else if (vm0 == LTR_MODEL_ARRHENIUS) { + tempDepType_ = 1; + } else if (vm0 == LTR_MODEL_NOTSET) { + throw CanteraError("SimpleTransport::initLiquid", + "Viscosity Model is not set in the input file"); + } else { + throw CanteraError("SimpleTransport::initLiquid", + "Viscosity Model is not handled by this object"); + } + + for (k = 0; k < m_nsp; k++) { + Cantera::LiquidTransportData <d = tr.LTData[k]; + LiquidTR_Model vm = ltd.model_viscosity; + if (vm != vm0) { + throw CanteraError(" SimpleTransport::initLiquid", + "different viscosity models"); + } + vector_fp &kentry = m_coeffVisc_Ns[k]; + kentry = ltd.viscCoeffs; + } + + /* + * Get the input thermal conductivities + */ + m_condSpecies.resize(m_nsp); + m_coeffLambda_Ns.clear(); + m_coeffLambda_Ns.resize(m_nsp); + LiquidTR_Model cm0 = ltd0.model_thermalCond; + if (cm0 != vm0) { + throw CanteraError("SimpleTransport::initLiquid", + "Conductivity model is not the same as the viscosity model"); + } + + for (k = 0; k < m_nsp; k++) { + Cantera::LiquidTransportData <d = tr.LTData[k]; + LiquidTR_Model cm = ltd.model_thermalCond; + if (cm != cm0) { + throw CanteraError(" SimpleTransport::initLiquid", + "different thermal conductivity models"); + } + vector_fp &kentry = m_coeffLambda_Ns[k]; + kentry = ltd.thermalCondCoeffs; + } + + /* + * Get the input species diffusivities + */ + useHydroRadius_ = false; + + m_diffSpecies.resize(m_nsp); + m_coeffDiff_Ns.clear(); + m_coeffDiff_Ns.resize(m_nsp); + LiquidTR_Model dm0 = ltd0.model_speciesDiffusivity; + if (dm0 != vm0) { + if (dm0 == LTR_MODEL_NOTSET) { + LiquidTR_Model rm0 = ltd0.model_hydroradius; + if (rm0 != vm0) { + throw CanteraError("SimpleTransport::initLiquid", + "hydroradius model is not the same as the viscosity model"); + } else { + useHydroRadius_ = true; + } + } + + for (k = 0; k < m_nsp; k++) { + Cantera::LiquidTransportData <d = tr.LTData[k]; + LiquidTR_Model dm = ltd.model_speciesDiffusivity; + if (dm == LTR_MODEL_NOTSET) { + LiquidTR_Model rm = ltd.model_hydroradius; + if (rm != vm0) { + throw CanteraError("SimpleTransport::initLiquid", + "hydroradius model is not the same as the viscosity model"); + } + if (rm != LTR_MODEL_CONSTANT) { + throw CanteraError("SimpleTransport::initLiquid", + "hydroradius model is not constant"); + } + vector_fp &kentry = m_coeffHydroRadius_Ns[k]; + kentry.push_back(ltd.hydroradius); + } else { + if (dm != dm0) { + throw CanteraError(" SimpleTransport::initLiquid", + "different thermal conductivity models"); + } + vector_fp &kentry = m_coeffDiff_Ns[k]; + kentry = ltd.speciesDiffusivityCoeffs; + } + } + + } + + + m_molefracs.resize(m_nsp); + m_concentrations.resize(m_nsp); + + m_chargeSpecies.resize(m_nsp); + for (k = 0; k < m_nsp; k++) { + m_chargeSpecies[k] = m_thermo->charge(k); + } + m_spwork.resize(m_nsp); + + // resize the internal gradient variables + m_Grad_X.resize(m_nDim * m_nsp, 0.0); + m_Grad_T.resize(m_nDim, 0.0); + m_Grad_P.resize(m_nDim, 0.0); + m_Grad_V.resize(m_nDim, 0.0); + + + + // set all flags to false + m_visc_mix_ok = false; + m_visc_temp_ok = false; + + m_cond_temp_ok = false; + m_cond_mix_ok = false; + + m_diff_temp_ok = false; + m_diff_mix_ok = false; + + return true; + } + + //================================================================================================ + // Returns the mixture viscosity of the solution + /* + * The viscosity is computed using the general mixture rules + * specified in the variable compositionDepType_. + * + * Solvent-only: + * \f[ + * \mu = \mu_0 + * \f] + * Mixture-average: + * \f[ + * \mu = \sum_k {\mu_k X_k} + * \f] + * + * Here \f$ \mu_k \f$ is the viscosity of pure species \e k. + * + * @see updateViscosity_T(); + */ + doublereal SimpleTransport::viscosity() { + + update_T(); + update_C(); + + if (m_visc_mix_ok) return m_viscmix; + + // update m_viscSpecies[] if necessary + if (!m_visc_temp_ok) { + updateViscosity_T(); + } + + if (compositionDepType_ == 0) { + m_viscmix = m_viscSpecies[0]; + } else if (compositionDepType_ == 1) { + m_viscmix = 0.0; + for (int k = 0; k < m_nsp; k++) { + m_viscmix += m_viscSpecies[k] * m_molefracs[k]; + } + } + m_visc_mix_ok = true; + return m_viscmix; + } + //================================================================================================ + void SimpleTransport::getSpeciesViscosities(doublereal* visc) { + update_T(); + if (!m_visc_temp_ok) { + updateViscosity_T(); + } + copy(m_viscSpecies.begin(), m_viscSpecies.end(), visc); + } + //================================================================================================ + void SimpleTransport::getBinaryDiffCoeffs(int ld, doublereal* d) { + int i, j; + double bdiff; + update_T(); + + // if necessary, evaluate the species diffusion coefficents + // from the polynomial fits + if (!m_diff_temp_ok) updateDiff_T(); + + for (i = 0; i < m_nsp; i++) { + for (j = 0; j < m_nsp; j++) { + bdiff = 0.5 * (m_diffSpecies[i] + m_diffSpecies[j]); + d[i*m_nsp+j] = bdiff; + } + } + } + //================================================================================================ + // Get the electrical Mobilities (m^2/V/s). + /* + * This function returns the mobilities. In some formulations + * this is equal to the normal mobility multiplied by faraday's constant. + * + * Frequently, but not always, the mobility is calculated from the + * diffusion coefficient using the Einstein relation + * + * \f[ + * \mu^e_k = \frac{F D_k}{R T} + * \f] + * + * @param mobil_e Returns the mobilities of + * the species in array \c mobil_e. The array must be + * dimensioned at least as large as the number of species. + */ + void SimpleTransport::getMobilities(doublereal* const mobil) { + int k; + getMixDiffCoeffs(DATA_PTR(m_spwork)); + doublereal c1 = ElectronCharge / (Boltzmann * m_temp); + for (k = 0; k < m_nsp; k++) { + mobil[k] = c1 * m_spwork[k]; + } + } + //================================================================================================ + // Get the fluid mobilities (s kmol/kg). + /* + * This function returns the fluid mobilities. Usually, you have + * to multiply Faraday's constant into the resulting expression + * to general a species flux expression. + * + * Frequently, but not always, the mobility is calculated from the + * diffusion coefficient using the Einstein relation + * + * \f[ + * \mu^f_k = \frac{D_k}{R T} + * \f] + * + * + * @param mobil_f Returns the mobilities of + * the species in array \c mobil. The array must be + * dimensioned at least as large as the number of species. + */ + void SimpleTransport::getFluidMobilities(doublereal* const mobil_f) { + int k; + getMixDiffCoeffs(DATA_PTR(m_spwork)); + doublereal c1 = 1.0 / (GasConstant * m_temp); + for (k = 0; k < m_nsp; k++) { + mobil_f[k] = c1 * m_spwork[k]; + } + } + //================================================================================================ + void SimpleTransport::set_Grad_V(const doublereal* const grad_V) { + doMigration_ = false; + for (int a = 0; a < m_nDim; a++) { + m_Grad_V[a] = grad_V[a]; + if (fabs(grad_V[a]) > 1.0E-13) doMigration_ = true; + } + } + //================================================================================================ + void SimpleTransport::set_Grad_T(const doublereal* const grad_T) { + for (int a = 0; a < m_nDim; a++) { + m_Grad_T[a] = grad_T[a]; + } + } + //================================================================================================ + void SimpleTransport::set_Grad_X(const doublereal* const grad_X) { + int itop = m_nDim * m_nsp; + for (int i = 0; i < itop; i++) { + m_Grad_X[i] = grad_X[i]; + } + } + //================================================================================================ + // Returns the mixture thermal conductivity of the solution + /* + * The thermal is computed using the general mixture rules + * specified in the variable compositionDepType_. + * + * Solvent-only: + * \f[ + * \lambda = \lambda_0 + * \f] + * Mixture-average: + * \f[ + * \lambda = \sum_k {\lambda_k X_k} + * \f] + * + * Here \f$ \lambda_k \f$ is the thermal conductivity of pure species \e k. + * + * @see updateCond_T(); + */ + doublereal SimpleTransport::thermalConductivity() { + update_T(); + update_C(); + if (!m_cond_temp_ok) { + updateCond_T(); + } + if (!m_cond_mix_ok) { + if (compositionDepType_ == 0) { + m_lambda = m_condSpecies[0]; + } else if (compositionDepType_ == 1) { + m_lambda = 0.0; + for (int k = 0; k < m_nsp; k++) { + m_lambda += m_condSpecies[k] * m_molefracs[k]; + } + } + m_cond_mix_ok = true; + } + return m_lambda; + } + //================================================================================================ + + /* + * Thermal diffusion is not considered in this mixture-averaged + * model. To include thermal diffusion, use transport manager + * MultiTransport instead. This methods fills out array dt with + * zeros. + */ + void SimpleTransport::getThermalDiffCoeffs(doublereal* const dt) { + for (int k = 0; k < m_nsp; k++) { + dt[k] = 0.0; + } + } +//================================================================================================ + /** + * @param ndim The number of spatial dimensions (1, 2, or 3). + * @param grad_T The temperature gradient (ignored in this model). + * @param ldx Leading dimension of the grad_X array. + * The diffusive mass flux of species \e k is computed from + * + * \f[ + * \vec{j}_k = -n M_k D_k \nabla X_k. + * \f] + */ + void SimpleTransport::getSpeciesFluxes(int ndim, + const doublereal* grad_T, + int ldx, const doublereal* grad_X, + int ldf, doublereal* fluxes) { + set_Grad_T(grad_T); + set_Grad_X(grad_X); + getSpeciesFluxesExt(ldf, fluxes); + } + //================================================================================================ + // Return the species diffusive mass fluxes wrt to + // the mass averaged velocity. + /* + * + * units = kg/m2/s + * + * Internally, gradients in the in mole fraction, temperature + * and electrostatic potential contribute to the diffusive flux + * + * + * The diffusive mass flux of species \e k is computed from the following + * formula + * + * \f[ + * j_k = - M_k z_k u^f_k F c_k \nabla \Psi - c M_k D_k \nabla X_k - Y_k V_c + * \f] + * + * where V_c is the correction velocity + * + * \f[ + * V_c = - \sum_j {M_k z_k u^f_k F c_k \nabla \Psi + c M_j D_j \nabla X_j} + * \f] + * + * @param ldf stride of the fluxes array. Must be equal to + * or greater than the number of species. + * @param fluxes Vector of calculated fluxes + */ + void SimpleTransport::getSpeciesFluxesExt(int ldf, doublereal* fluxes) { + int n, k; + AssertThrow(ldf >= m_nsp ,"SimpleTransport::getSpeciesFluxesExt: Stride must be greater than m_nsp"); + update_T(); + update_C(); + + getMixDiffCoeffs(DATA_PTR(m_spwork)); + + const array_fp& mw = m_thermo->molecularWeights(); + const doublereal* y = m_thermo->massFractions(); + doublereal conc = m_thermo->molarDensity(); + // Unroll wrt ndim + + vector_fp sum(m_nDim, 0.0); + + if (doMigration_) { + for (n = 0; n < m_nDim; n++) { + for (k = 0; k < m_nsp; k++) { + fluxes[n*ldf + k] = -conc * mw[k] * m_spwork[k] * m_Grad_X[n*m_nsp + k]; + sum[n] += fluxes[n*ldf + k]; + } + } + } else { + double FRT = ElectronCharge / (Boltzmann * m_temp); + for (n = 0; n < m_nDim; n++) { + for (k = 0; k < m_nsp; k++) { + fluxes[n*ldf + k] = -conc * mw[k] * m_spwork[k] * + ( m_Grad_X[n*m_nsp + k] + FRT * m_molefracs[k] * m_chargeSpecies[k] * m_Grad_V[n*m_nsp + k]); + sum[n] += fluxes[n*ldf + k]; + } + } + } + // add correction flux to enforce sum to zero + for (n = 0; n < m_nDim; n++) { + for (k = 0; k < m_nsp; k++) { + fluxes[n*ldf + k] -= y[k]*sum[n]; + } + } + } + //================================================================================================ + // Mixture-averaged diffusion coefficients [m^2/s]. + /* + * Returns the simple diffusion coefficients input into the model. Nothing fancy here. + */ + void SimpleTransport::getMixDiffCoeffs(doublereal* const d) { + update_T(); + update_C(); + // update the binary diffusion coefficients if necessary + if (!m_diff_temp_ok) { + updateDiff_T(); + } + for (int k = 0; k < m_nsp; k++) { + d[k] = m_diffSpecies[k]; + } + } + //================================================================================================ + + // Handles the effects of changes in the mixture concentration + /* + * This is called for every interface call to check whether + * the concentrations have changed. Concentrations change + * whenever the pressure or the mole fraction has changed. + * If it has changed, the recalculations should be done. + * + * Note this should be a lightweight function since it's + * part of all of the interfaces. + * + * @internal + */ + bool SimpleTransport::update_C() { + // If the pressure has changed then the concentrations + // have changed. + doublereal pres = m_thermo->pressure(); + bool qReturn = true; + if (pres != m_press) { + qReturn = false; + m_press = pres; + } + int iStateNew = m_thermo->stateMFNumber(); + if (iStateNew != m_iStateMF) { + qReturn = false; + m_thermo->getMoleFractions(DATA_PTR(m_molefracs)); + m_thermo->getConcentrations(DATA_PTR(m_concentrations)); + concTot_ = 0.0; + for (int k = 0; k < m_nsp; k++) { + m_molefracs[k] = fmaxx(0.0, m_molefracs[k]); + concTot_ += m_concentrations[k]; + } + dens_ = m_thermo->density(); + meanMolecularWeight_ = m_thermo->meanMolecularWeight(); + } + if (qReturn) { + return false; + } + + + // Mixture stuff needs to be evaluated + m_visc_mix_ok = false; + m_diff_mix_ok = false; + m_cond_mix_ok = false; + + return true; + } + + //================================================================================================ + /** + * Update the temperature-dependent parts of the mixture-averaged + * thermal conductivity. + */ + void SimpleTransport::updateCond_T() { + int k; + if (tempDepType_ == 0) { + for (k = 0; k < m_nsp; k++) { + Coeff_T_ &coeff = m_coeffLambda_Ns[k]; + m_condSpecies[k] = coeff[0]; + } + } else if (tempDepType_ == 1) { + for (k = 0; k < m_nsp; k++) { + Coeff_T_ &coeff = m_coeffLambda_Ns[k]; + m_condSpecies[k] = coeff[0] * pow(m_temp,coeff[1]) * exp(-coeff[2]/m_temp); + } + } + m_cond_temp_ok = true; + m_cond_mix_ok = false; + } + //================================================================================================ + /** + * Update the species diffusion coefficients. + */ + void SimpleTransport::updateDiff_T() { + int k; + if (useHydroRadius_) { + if (tempDepType_ == 0) { + for (k = 0; k < m_nsp; k++) { + Coeff_T_ &coeff = m_coeffDiff_Ns[k]; + m_diffSpecies[k] = coeff[0]; + } + } else if (tempDepType_ == 1) { + for (k = 0; k < m_nsp; k++) { + Coeff_T_ &coeff = m_coeffDiff_Ns[k]; + m_diffSpecies[k] = coeff[0] * pow(m_temp,coeff[1]) * exp(-coeff[2]/m_temp); + } + } + } else { + double visc = viscosity(); + double RT = GasConstant * m_temp; + for (k = 0; k < m_nsp; k++) { + Coeff_T_ &coeff = m_coeffHydroRadius_Ns[k]; + double rad = coeff[0]; + m_diffSpecies[k] = RT / (6.0 * Pi * visc * rad); + } + } + m_diff_temp_ok = true; + m_diff_mix_ok = false; + } + //================================================================================================ + /** + * Update the pure-species viscosities. + */ + void SimpleTransport::updateViscosities_C() { + + } + //================================================================================================ + /** + * Update the temperature-dependent viscosity terms. + * Updates the array of pure species viscosities, and the + * weighting functions in the viscosity mixture rule. + * The flag m_visc_ok is set to true. + */ + void SimpleTransport::updateViscosity_T() { + int k; + if (tempDepType_ == 0) { + for (k = 0; k < m_nsp; k++) { + Coeff_T_ &coeff = m_coeffVisc_Ns[k]; + m_viscSpecies[k] = coeff[0]; + } + } else if (tempDepType_ == 1) { + for (k = 0; k < m_nsp; k++) { + Coeff_T_ &coeff = m_coeffVisc_Ns[k]; + m_viscSpecies[k] = coeff[0] * pow(m_temp,coeff[1]) * exp(-coeff[2]/m_temp); + } + } + m_visc_temp_ok = true; + m_visc_mix_ok = false; + } + //================================================================================================= + bool SimpleTransport::update_T() + { + doublereal t = m_thermo->temperature(); + if (t == m_temp) return false; + if (t < 0.0) { + throw CanteraError("SimpleTransport::update_T", + "negative temperature "+fp2str(t)); + } + + // Compute various functions of temperature + m_temp = t; + + // temperature has changed, so polynomial temperature + // interpolations will need to be reevaluated. + // Set all of these flags to false + m_visc_mix_ok = false; + m_visc_temp_ok = false; + + m_cond_temp_ok = true; + m_cond_mix_ok = false; + + m_diff_mix_ok = false; + m_diff_temp_ok = false; + + return true; + } + //================================================================================================ + /** + * Throw an exception if this method is invoked. + * This probably indicates something is not yet implemented. + */ + doublereal SimpleTransport::err(std::string msg) const { + throw CanteraError("SimpleTransport Class", + "\n\n\n**** Method "+ msg +" not implemented in model " + + int2str(model()) + " ****\n" + "(Did you forget to specify a transport model?)\n\n\n"); + + return 0.0; + } + //================================================================================================ + +} +//================================================================================================ diff --git a/Cantera/src/transport/SimpleTransport.h b/Cantera/src/transport/SimpleTransport.h new file mode 100644 index 000000000..35eb9560f --- /dev/null +++ b/Cantera/src/transport/SimpleTransport.h @@ -0,0 +1,718 @@ +/** + * + * @file SimpleTransport.h + * Header file defining class SimpleTransport + */ +/* + * $Revision: 1.9 $ + * $Date: 2009/03/27 18:24:39 $ + */ + +#ifndef CT_SIMPLETRAN_H +#define CT_SIMPLETRAN_H + + + +// STL includes +#include +#include +#include +#include +#include + +using namespace std; + +// Cantera includes +#include "TransportBase.h" +#include "DenseMatrix.h" +#include "TransportParams.h" +#include "LiquidTransportParams.h" + +namespace Cantera { + + + + + class LiquidTransportParams; + + + //! Class LiquidTransport implements mixture-averaged transport + //! properties for liquid phases. + /*! + * The model is based on that + * described by Newman, Electrochemical Systems + * + * The velocity of species i may be described by the + * following equation p. 297 (12.1) + * + * \f[ + * c_i \nabla \mu_i = R T \sum_j \frac{c_i c_j}{c_T D_{ij}} + * (\mathbf{v}_j - \mathbf{v}_i) + * \f] + * + * This as written is degenerate by 1 dof. + * + * To fix this we must add in the definition of the mass averaged + * velocity of the solution. We will call the simple bold-faced + * \f$\mathbf{v} \f$ + * symbol the mass-averaged velocity. Then, the relation + * between \f$\mathbf{v}\f$ and the individual species velocities is + * \f$\mathbf{v}_i\f$ + * + * \f[ + * \rho_i \mathbf{v}_i = \rho_i \mathbf{v} + \mathbf{j}_i + * \f] + * where \f$\mathbf{j}_i\f$ are the diffusional fluxes of species i + * with respect to the mass averaged velocity and + * + * \f[ + * \sum_i \mathbf{j}_i = 0 + * \f] + * + * and + * + * \f[ + * \sum_i \rho_i \mathbf{v}_i = \rho \mathbf{v} + * \f] + * + * Using these definitions, we can write + * + * \f[ + * \mathbf{v}_i = \mathbf{v} + \frac{\mathbf{j}_i}{\rho_i} + * \f] + * + * + * \f[ + * c_i \nabla \mu_i = R T \sum_j \frac{c_i c_j}{c_T D_{ij}} + * (\frac{\mathbf{j}_j}{\rho_j} - \frac{\mathbf{j}_i}{\rho_i}) + * = R T \sum_j \frac{1}{D_{ij}} + * (\frac{x_i \mathbf{j}_j}{M_j} - \frac{x_j \mathbf{j}_i}{M_i}) + * \f] + * + * The equations that we actually solve are + * + * \f[ + * c_i \nabla \mu_i = + * = R T \sum_j \frac{1}{D_{ij}} + * (\frac{x_i \mathbf{j}_j}{M_j} - \frac{x_j \mathbf{j}_i}{M_i}) + * \f] + * and we replace the 0th equation with the following: + * + * \f[ + * \sum_i \mathbf{j}_i = 0 + * \f] + * + * When there are charged species, we replace the rhs with the + * gradient of the electrochemical potential to obtain the + * modified equation + * + * \f[ + * c_i \nabla \mu_i + c_i F z_i \nabla \Phi + * = R T \sum_j \frac{1}{D_{ij}} + * (\frac{x_i \mathbf{j}_j}{M_j} - \frac{x_j \mathbf{j}_i}{M_i}) + * \f] + * + * With this formulation we may solve for the diffusion velocities, + * without having to worry about what the mass averaged velocity + * is. + * + *

Viscosity Calculation

+ * + * The viscosity calculation may be broken down into two parts. + * In the first part, the viscosity of the pure species are calculated + * In the second part, a mixing rule is applied, based on the + * Wilkes correlation, to yield the mixture viscosity. + * + * + * + */ + class SimpleTransport : public Transport { + public: + + typedef vector_fp Coeff_T_; + + + //! Default constructor. + /*! + * This requires call to initLiquid(LiquidTransportParams& tr) + * after filling LiquidTransportParams to complete instantiation. + * The filling of LiquidTransportParams is currently carried out + * in the TransportFactory class, but might be moved at some point. + * + * @param thermo ThermoPhase object holding species information. + * @param ndim Number of spatial dimensions. + */ + SimpleTransport(thermo_t* thermo = 0, int ndim = 1); + + //!Copy Constructor for the %LiquidThermo object. + /*! + * @param right %LiquidTransport to be copied + */ + SimpleTransport(const SimpleTransport &right); + + //! Assignment operator + /*! + * This is NOT a virtual function. + * + * @param right Reference to %LiquidTransport object to be copied + * into the current one. + */ + SimpleTransport& operator=(const SimpleTransport& right); + + //! Duplication routine for objects which inherit from + //! %Transport + /*! + * This virtual routine can be used to duplicate %Transport objects + * inherited from %Transport even if the application only has + * a pointer to %Transport to work with. + * + * These routines are basically wrappers around the derived copy + * constructor. + */ + virtual Transport *duplMyselfAsTransport() const; + + + //! virtual destructor + virtual ~SimpleTransport() {} + + //! Initialize the transport object + /*! + * Here we change all of the internal dimensions to be sufficient. + * We get the object ready to do property evaluations. + * + * @param tr Transport parameters for all of the species + * in the phase. + */ + virtual bool initLiquid(LiquidTransportParams& tr); + + friend class TransportFactory; + + + //! Return the model id for this transport parameterization + virtual int model() const { + return cSimpleTransport; + } + + //! overloaded base class methods + + //! Returns the mixture viscosity of the solution + /*! + * The viscosity is computed using the general mixture rules + * specified in the variable compositionDepType_. + * + * Solvent-only: + * \f[ + * \mu = \mu_0 + * \f] + * Mixture-average: + * \f[ + * \mu = \sum_k {\mu_k X_k} + * \f] + * + * Here \f$ \mu_k \f$ is the viscosity of pure species \e k. + * + * @see updateViscosity_T(); + */ + virtual doublereal viscosity(); + + //! Returns the pure species viscosities + /*! + * The pure species viscosities are to be given in an Arrhenius + * form in accordance with activated-jump-process dominated transport. + */ + virtual void getSpeciesViscosities(doublereal* const visc); + + //! Returns the binary diffusion coefficients + /*! + * @param ld + * @param d + */ + virtual void getBinaryDiffCoeffs(const int ld, doublereal* const d); + + //! Get the Mixture diffusion coefficients + /*! + * @param d vector of mixture diffusion coefficients + * units = m2 s-1. length = number of species + */ + virtual void getMixDiffCoeffs(doublereal* const d); + + + //! Return the thermal diffusion coefficients + /*! + * These are all zero for this simple implementaion + * + * @param dt thermal diffusion coefficients + */ + virtual void getThermalDiffCoeffs(doublereal* const dt); + + + //! Returns the mixture thermal conductivity of the solution + /*! + * The thermal is computed using the general mixture rules + * specified in the variable compositionDepType_. + * + * Controlling update boolean = m_condmix_ok + * + * Units are in W/m/K or equivalently kg m / s3 / K + * + * Solvent-only: + * \f[ + * \lambda = \lambda_0 + + * \f] + * Mixture-average: + * \f[ + * \lambda = \sum_k {\lambda_k X_k} + * \f] + * + * Here \f$ \lambda_k \f$ is the thermal conductivity of pure species \e k. + * + * @see updateCond_T(); + */ + + virtual doublereal thermalConductivity(); + + //! Get the electrical Mobilities (m^2/V/s). + /*! + * This function returns the mobilities. In some formulations + * this is equal to the normal mobility multiplied by faraday's constant. + * + * Frequently, but not always, the mobility is calculated from the + * diffusion coefficient using the Einstein relation + * + * \f[ + * \mu^e_k = \frac{F D_k}{R T} + * \f] + * + * @param mobil_e Returns the mobilities of + * the species in array \c mobil_e. The array must be + * dimensioned at least as large as the number of species. + */ + virtual void getMobilities(doublereal* const mobil_e); + + //! Get the fluid mobilities (s kmol/kg). + /*! + * This function returns the fluid mobilities. Usually, you have + * to multiply Faraday's constant into the resulting expression + * to general a species flux expression. + * + * Frequently, but not always, the mobility is calculated from the + * diffusion coefficient using the Einstein relation + * + * \f[ + * \mu^f_k = \frac{D_k}{R T} + * \f] + * + * + * @param mobil_f Returns the mobilities of + * the species in array \c mobil. The array must be + * dimensioned at least as large as the number of species. + */ + virtual void getFluidMobilities(doublereal* const mobil_f); + + //! Specify the valpdaue of the gradient of the voltage + /*! + * + * @param grad_V Gradient of the voltage (length num dimensions); + */ + virtual void set_Grad_V(const doublereal* const grad_V); + + //! Specify the value of the gradient of the temperature + /*! + * + * @param grad_V Gradient of the temperature (length num dimensions); + */ + virtual void set_Grad_T(const doublereal* const grad_T); + + //! Specify the value of the gradient of the MoleFractions + /*! + * + * @param grad_X Gradient of the mole fractions(length nsp * num dimensions); + */ + virtual void set_Grad_X(const doublereal* const grad_X); + + + /** + * @param ndim The number of spatial dimensions (1, 2, or 3). + * @param grad_T The temperature gradient (ignored in this model). + * @param ldx Leading dimension of the grad_X array. + * The diffusive mass flux of species \e k is computed from + * + * + */ + virtual void getSpeciesFluxes(int ndim, + const doublereal* grad_T, + int ldx, const doublereal* grad_X, + int ldf, doublereal* fluxes); + + //! Return the species diffusive mass fluxes wrt to + //! the mass averaged velocity, + /*! + * + * units = kg/m2/s + * + * Internally, gradients in the in mole fraction, temperature + * and electrostatic potential contribute to the diffusive flux + * + * + * The diffusive mass flux of species \e k is computed from the following + * formula + * + * \f[ + * j_k = - \rho M_k D_k \nabla X_k - Y_k V_c + * \f] + * + * where V_c is the correction velocity + * + * \f[ + * V_c = - \sum_j {\rho M_j D_j \nabla X_j} + * \f] + * + * @param ldf stride of the fluxes array. Must be equal to + * or greater than the number of species. + * @param fluxes Vector of calculated fluxes + */ + virtual void getSpeciesFluxesExt(int ldf, doublereal* fluxes); + + protected: + + //! Handles the effects of changes in the Temperature, internally + //! within the object. + /*! + * This is called whenever a transport property is requested. + * The first task is to check whether the temperature has changed + * since the last call to update_T(). + * If it hasn't then an immediate return is carried out. + * + * @internal + * + * @return Returns true if the temperature has changed, and false otherwise + */ + virtual bool update_T(); + + //! Handles the effects of changes in the mixture concentration + /*! + * This is called for every interface call to check whether + * the concentrations have changed. Concentrations change + * whenever the pressure or the mole fraction has changed. + * If it has changed, the recalculations should be done. + * + * Note this should be a lightweight function since it's + * part of all of the interfaces. + * + * @internal + */ + virtual bool update_C(); + + //! Update the temperature-dependent viscosity terms. + //! Updates the array of pure species viscosities, and the + //! weighting functions in the viscosity mixture rule. + /*! + * The flag m_visc_temp_ok is set to true. + */ + void updateViscosity_T(); + + //! Update the temperature-dependent parts of the mixture-averaged + //! thermal conductivity. + void updateCond_T(); + + //! Update the concentration parts of the viscosities + /*! + * Internal routine is run whenever the update_boolean + * is false. This routine will calculate + * internal values for the species viscosities. + * + * @internal + */ + void updateViscosities_C(); + + //! Update the binary diffusion coefficients wrt T. + /*! + * These are evaluated + * from the polynomial fits at unit pressure (1 Pa). + */ + void updateDiff_T(); + + + private: + + //! Number of species in the mixture + int m_nsp; + + //! Temperature dependence type + /*! + * The following coefficients are allowed to have simple + * temperature dependencies: + * mixture viscosity + * mixture thermal conductivity + * diffusitivy + * + * Types of temperature dependencies: + * 0 - Independent of temperature (only one implemented so far) + * 1 - extended arrhenius form + * 2 - power law form + */ + int tempDepType_; + + //! Composition dependence of the transport properties + /*! + * The following coefficients are allowed to have simple + * composition dependencies + * mixture viscosity + * mixture thermal conductivity + * + * + * Types of composition dependencies + * 0 - Solvent values (i.e., species 0) contributes only + * 1 - linear combination of mole fractions; + */ + int compositionDepType_; + + bool useHydroRadius_; + + //! Boolean indicating whether electro-migration term should be + //! added + /*! + * + */ + bool doMigration_; + + //! Minimum temperature applicable to the transport property eval + doublereal m_tmin; + + //! Maximum temperature applicable to the transport property evaluator + doublereal m_tmax; + + //! Local Copy of the molecular weights of the species + /*! + * Length is Equal to the number of species in the mechanism. + */ + vector_fp m_mw; + + //! Pure species viscosities in Arrhenius temperature-dependent form. + std::vector m_coeffVisc_Ns; + + //! Pure species thermal conductivities in Arrhenius temperature-dependent form. + /*! + * + */ + std::vector m_coeffLambda_Ns; + + + //! Pure species viscosities in Arrhenius temperature-dependent form. + std::vector m_coeffDiff_Ns; + + + std::vector m_coeffHydroRadius_Ns; + + + //! Internal value of the gradient of the mole fraction vector + /*! + * Note, this is the only gradient value that can and perhaps + * should reflect the true state of the mole fractions in the + * application solution vector. In other words no cropping or + * massaging of the values to make sure they are above zero + * should occur. - developing .... + * + * m_nsp is the number of species in the fluid + * k is the species index + * n is the dimensional index (x, y, or z). It has a length + * equal to m_nDim + * + * m_Grad_X[n*m_nsp + k] + */ + vector_fp m_Grad_X; + + //! Internal value of the gradient of the Temperature vector + /*! + * Generally, if a transport property needs this + * in its evaluation it will look to this place + * to get it. + * + * No internal property is precalculated based on gradients. + * Gradients are assumed to be freshly updated before + * every property call. + */ + vector_fp m_Grad_T; + + //! Internal value of the gradient of the Pressure vector + /*! + * Generally, if a transport property needs this + * in its evaluation it will look to this place + * to get it. + * + * No internal property is precalculated based on gradients. + * Gradients are assumed to be freshly updated before + * every property call. + */ + vector_fp m_Grad_P; + + //! Internal value of the gradient of the Electric Voltage + /*! + * Generally, if a transport property needs this + * in its evaluation it will look to this place + * to get it. + * + * No internal property is precalculated based on gradients. + * Gradients are assumed to be freshly updated before + * every property call. + */ + vector_fp m_Grad_V; + + + // property values + + + + //! Vector of Species Diffusivities + /*! + * Depends on the temperature. We have set the pressure dependence + * to zero for this liquid phase constituitve model + * + * units m2/s + */ + vector_fp m_diffSpecies; + + //! Species viscosities + /*! + * Viscosity of the species + * Length = number of species + * + * Depends on the temperature. We have set the pressure dependence + * to zero for this model + * + * controlling update boolean -> m_visc_temp_ok + */ + vector_fp m_viscSpecies; + + //! Internal value of the species individual thermal conductivities + /*! + * Then a mixture rule is applied to get the solution conductivities + * + * Depends on the temperature and perhaps pressure, but + * not the species concentrations + * + * controlling update boolean -> m_cond_temp_ok + */ + vector_fp m_condSpecies; + + //! State of the mole fraction vector. + int m_iStateMF; + + //! Local copy of the mole fractions of the species in the phase + /*! + * The mole fractions here are assumed to be bounded by 0.0 and 1.0 + * and they are assumed to add up to one exactly. This mole + * fraction vector comes from the ThermoPhase object. Derivative + * quantities from this are referred to as bounded. + * + * Update info? + * length = m_nsp + */ + vector_fp m_molefracs; + + + //! Local copy of the concentrations of the species in the phase + /*! + * The concentrations are consistent with the m_molefracs + * vector which is bounded and sums to one. + * + * Update info? + * length = m_nsp + */ + vector_fp m_concentrations; + + //! Local copy of the total concentration. + /*! + * This is consistent with the m_concentrations[] and + * m_molefracs[] vector. + */ + doublereal concTot_; + + //! Mean molecular weight + doublereal meanMolecularWeight_; + + //! Density + doublereal dens_; + + //! Local copy of the charge of each species + /*! + * Contains the charge of each species (length m_nsp) + */ + vector_fp m_chargeSpecies; + + //! Current Temperature -> locally storred + /*! + * This is used to test whether new temperature computations + * should be performed. + */ + doublereal m_temp; + + + //! Current value of the pressure + doublereal m_press; + + + //! Saved value of the mixture thermal conductivity + doublereal m_lambda; + + //! Saved value of the mixture viscosity + doublereal m_viscmix; + + //! work space + /*! + * Length is equal to m_nsp + */ + vector_fp m_spwork; + + + + private: + //! Boolean indicating that the top-level mixture viscosity is current + /*! + * This is turned false for every change in T, P, or C. + */ + bool m_visc_mix_ok; + + //! Boolean indicating that weight factors wrt viscosity is current + bool m_visc_temp_ok; + + //! Boolean indicating that mixture diffusion coeffs are current + bool m_diff_mix_ok; + + //! Boolean indicating that binary diffusion coeffs are current + bool m_diff_temp_ok; + + //! Flag to indicate that the pure species conductivities + //! are current wrt the temperature + bool m_cond_temp_ok; + + //! Boolean indicating that mixture conductivity is current + bool m_cond_mix_ok; + + + //! Number of dimensions + /*! + * Either 1, 2, or 3 + */ + int m_nDim; + + private: + + //! Throw an exception if this method is invoked. + /*! + * This probably indicates something is not yet implemented. + * + * @pram msg Indicates the member function which is not implemented + */ + doublereal err(std::string msg) const; + + }; +} +#endif + + + + + + diff --git a/Cantera/src/transport/SolidTransport.h b/Cantera/src/transport/SolidTransport.h index 027ba9816..76e18f943 100644 --- a/Cantera/src/transport/SolidTransport.h +++ b/Cantera/src/transport/SolidTransport.h @@ -46,7 +46,7 @@ namespace Cantera { public: virtual ~SolidTransport() {} - virtual int model() { return cSolidTransport; } + virtual int model() const { return cSolidTransport; } virtual doublereal thermalConductivity(); virtual void getMixDiffCoeffs(doublereal* const d); diff --git a/Cantera/src/transport/TransportBase.cpp b/Cantera/src/transport/TransportBase.cpp index 532530690..75da215b7 100644 --- a/Cantera/src/transport/TransportBase.cpp +++ b/Cantera/src/transport/TransportBase.cpp @@ -9,10 +9,12 @@ #include "ThermoPhase.h" #include "LiquidTransport.h" +#include "ctexceptions.h" #include "utilities.h" #include "LiquidTransportParams.h" #include "TransportFactory.h" +#include "stringUtils.h" #include "ctlapack.h" @@ -31,8 +33,6 @@ namespace Cantera { //////////////////// class LiquidTransport methods ////////////// - - Transport::Transport(thermo_t* thermo, int ndim) : m_thermo(thermo), m_ready(false), @@ -109,5 +109,37 @@ namespace Cantera { { err("setParameters"); } + + + void Transport::setThermo(thermo_t& thermo) { + if (!ready()) { + m_thermo = &thermo; + m_nmin = m_thermo->nSpecies(); + } + else + throw CanteraError("Transport::setThermo", + "the phase object cannot be changed after " + "the transport manager has been constructed."); + } + + + doublereal Transport::err(std::string msg) const { + + throw CanteraError("Transport Base Class", + "\n\n\n**** Method "+ msg +" not implemented in model " + + int2str(model()) + " ****\n" + "(Did you forget to specify a transport model?)\n\n\n"); + + return 0.0; + } + + void Transport::finalize() { + if (!ready()) + m_ready = true; + else + throw CanteraError("Transport::finalize", + "finalize has already been called."); + } + } diff --git a/Cantera/src/transport/TransportBase.h b/Cantera/src/transport/TransportBase.h old mode 100755 new mode 100644 index 0ac33e1a7..864469f8c --- a/Cantera/src/transport/TransportBase.h +++ b/Cantera/src/transport/TransportBase.h @@ -28,6 +28,8 @@ namespace Cantera { class TransportParams; + class GasTransportParams; + class LiquidTransportParams; const int CK_Mode = 10; @@ -43,6 +45,7 @@ namespace Cantera { const int cFtnTransport = 600; const int cLiquidTransport = 700; const int cAqueousTransport = 750; + const int cSimpleTransport = 770; const int cRadiativeTransport = 800; const int cWaterTransport = 721; @@ -107,7 +110,7 @@ namespace Cantera { * virtual method returns an integer flag that identifies the * transport model implemented. The base class returns 0. */ - virtual int model() {return 0;} + virtual int model() const {return 0;} /** * Phase object. Every transport manager is designed to compute @@ -157,6 +160,14 @@ namespace Cantera { virtual doublereal viscosity() { return err("viscosity"); } + //! Returns the pure species viscosities + /*! + * The units are Pa-s and the length is the number of species + * + * @param visc Vector of viscosities + */ + virtual void getSpeciesViscosities(doublereal* const visc) + { err("getSpeciesViscosities"); } /** * The bulk viscosity in Pa-s. The bulk viscosity is only @@ -167,9 +178,11 @@ namespace Cantera { virtual doublereal bulkViscosity() { return err("bulkViscosity"); } - - /** - * The thermal conductivity in W/m/K. + //! Returns the mixture thermal conductivity in W/m/K. + /*! + * Units are in W / m K or equivalently kg m / s3 K + * + * @return returns thermal conductivity in W/m/K. */ virtual doublereal thermalConductivity() { return err("thermalConductivity"); } @@ -180,14 +193,48 @@ namespace Cantera { virtual doublereal electricalConductivity() { return err("electricalConductivity"); } - /** - * Electrical mobilities (m^2/V/s). Returns the mobilities of - * the species in array \c mobil. The array must be - * dimensioned at least as large as the number of species. + + //! Get the Electrical mobilities (m^2/V/s). + /*! + * This function returns the mobilities. In some formulations + * this is equal to the normal mobility multiplied by faraday's constant. + * + * Frequently, but not always, the mobility is calculated from the + * diffusion coefficient using the Einstein relation + * + * \f[ + * \mu^e_k = \frac{F D_k}{R T} + * \f] + * + * + * @param mobil_e Returns the mobilities of + * the species in array \c mobil_e. The array must be + * dimensioned at least as large as the number of species. */ - virtual void getMobilities(doublereal* const mobil) + virtual void getMobilities(doublereal* const mobil_e) { err("getMobilities"); } + //! Get the fluid mobilities (s kmol/kg). + /*! + * This function returns the fluid mobilities. Usually, you have + * to multiply Faraday's constant into the resulting expression + * to general a species flux expression. + * + * Frequently, but not always, the mobility is calculated from the + * diffusion coefficient using the Einstein relation + * + * \f[ + * \mu^f_k = \frac{D_k}{R T} + * \f] + * + * + * @param mobil_f Returns the mobilities of + * the species in array \c mobil. The array must be + * dimensioned at least as large as the number of species. + */ + virtual void getFluidMobilities(doublereal* const mobil_f) + { err("getFluidMobilities"); } + //@} @@ -221,6 +268,39 @@ namespace Cantera { err("getSpeciesFluxes"); } + //! Get the species diffusive mass fluxes wrt to + //! the mass averaged velocity, + //! given the gradients in mole fraction, temperature + //! and electrostatic potential. + /*! + * Units for the returned fluxes are kg m-2 s-1. + * + * @param ndim Number of dimensions in the flux expressions + * @param grad_T Gradient of the temperature + * (length = ndim) + * @param ldx Leading dimension of the grad_X array + * (usually equal to m_nsp but not always) + * @param grad_X Gradients of the mole fraction + * Flat vector with the m_nsp in the inner loop. + * length = ldx * ndim + * @param ldf Leading dimension of the fluxes array + * (usually equal to m_nsp but not always) + * @param grad_Phi Gradients of the electrostatic potential + * (length = ndim) + * @param fluxes Output of the diffusive mass fluxes + * Flat vector with the m_nsp in the inner loop. + * length = ldx * ndim + */ + virtual void getSpeciesFluxesES(int ndim, + const doublereal* grad_T, + int ldx, + const doublereal* grad_X, + int ldf, + const doublereal* grad_Phi, + doublereal* fluxes) { + getSpeciesFluxes( ndim, grad_T, ldx, grad_X, ldf, fluxes ); + } + /** * Get the molar fluxes [kmol/m^2/s], given the thermodynamic * state at two nearby points. @@ -238,9 +318,9 @@ namespace Cantera { * Get the mass fluxes [kg/m^2/s], given the thermodynamic * state at two nearby points. * @param state1 Array of temperature, density, and mass - * fractions for state 1. + * fractions for state 1. * @param state2 Array of temperature, density, and mass - * fractions for state 2. + * fractions for state 2. * @param delta Distance from state 1 to state 2 (m). */ virtual void getMassFluxes(const doublereal* state1, @@ -263,9 +343,13 @@ namespace Cantera { { err("getThermalDiffCoeffs"); } - /** - * Binary diffusion coefficients [m^2/s]. - */ + //! Returns the matrix of binary diffusion coefficients [m^2/s]. + /*! + * @param ld Inner stride for writing the two dimension diffusion + * coefficients into a one dimensional vector + * @param d Diffusion coefficient matrix (must be at least m_k * m_k + * in length. + */ virtual void getBinaryDiffCoeffs(const int ld, doublereal* const d) { err("getBinaryDiffCoeffs"); } @@ -313,8 +397,22 @@ namespace Cantera { /** * Called by TransportFactory to set parameters. */ - virtual bool init(TransportParams& tr) - { err("init"); return false; } + //virtual bool init(TransportParams& tr) + //{ err("init"); return false; } + + /** + * Called by TransportFactory to set parameters. + */ + virtual bool initGas( GasTransportParams& tr ) + { err("initGas"); return false; } + + /** + * Called by TransportFactory to set parameters. + */ + virtual bool initLiquid( LiquidTransportParams& tr ) + { err("initLiquid"); return false; } + + /** diff --git a/Cantera/src/transport/TransportFactory.cpp b/Cantera/src/transport/TransportFactory.cpp old mode 100755 new mode 100644 index fbdfce92f..3f409b8db --- a/Cantera/src/transport/TransportFactory.cpp +++ b/Cantera/src/transport/TransportFactory.cpp @@ -20,6 +20,7 @@ #include "MixTransport.h" #include "SolidTransport.h" #include "DustyGasTransport.h" +#include "SimpleTransport.h" #ifdef WITH_IDEAL_SOLUTIONS #include "LiquidTransport.h" @@ -95,31 +96,21 @@ namespace Cantera { #endif }; - //////////////////// class Transport methods ///////////////////// - void Transport::setThermo(thermo_t& thermo) { - if (!ready()) { - m_thermo = &thermo; - m_nmin = m_thermo->nSpecies(); - } - else - throw CanteraError("Transport::setThermo", - "the phase object cannot be changed after " - "the transport manager has been constructed."); - } - - void Transport::finalize() { - if (!ready()) - m_ready = true; - else - throw CanteraError("Transport::finalize", - "finalize has already been called."); - } - - doublereal Transport::err(string msg) const { - throw NotImplemented(msg); - //return 0.0; - } + /** + * getArrhenius() parses the xml element called Arrhenius. + * The Arrhenius expression is + * \f[ k = A T^(b) exp (-E_a / RT). \f] + */ + static void getArrhenius(const XML_Node& node, + doublereal& A, doublereal& b, doublereal& E) { + /* parse the children for the A, b, and E conponents. + */ + A = getFloat(node, "A", "toSI"); + b = getFloat(node, "b"); + E = getFloat(node, "E", "actEnergy"); + E /= GasConstant; + } //////////////////// class TransportFactory methods ////////////// @@ -148,7 +139,7 @@ namespace Cantera { * @note This method is not used currently. */ void TransportFactory::getBinDiffCorrection(doublereal t, - const TransportParams& tr, int k, int j, doublereal xk, doublereal xj, + const GasTransportParams& tr, int k, int j, doublereal xk, doublereal xj, doublereal& fkj, doublereal& fjk) { doublereal w1, w2, wsum, sig1, sig2, sig12, sigratio, sigratio2, @@ -220,7 +211,7 @@ namespace Cantera { * correction, see Dixon-Lewis, Proc. Royal Society (1968). */ void TransportFactory::makePolarCorrections(int i, int j, - const TransportParams& tr, doublereal& f_eps, doublereal& f_sigma) { + const GasTransportParams& tr, doublereal& f_eps, doublereal& f_sigma) { // no correction if both are nonpolar, or both are polar if (tr.polar[i] == tr.polar[j]) { @@ -262,6 +253,9 @@ namespace Cantera { m_models["DustyGas"] = cDustyGasTransport; m_models["CK_Multi"] = CK_Multicomponent; m_models["CK_Mix"] = CK_MixtureAveraged; + m_models["Liquid"] = cLiquidTransport; + m_models["Aqueous"] = cAqueousTransport; + m_models["Simple"] = cSimpleTransport; m_models["User"] = cUserTransport; m_models["None"] = None; //m_models["Radiative"] = cRadiative; @@ -301,7 +295,7 @@ namespace Cantera { * make one of several transport models, and return a base class * pointer to it. */ - Transport* TransportFactory::newTransport(string transportModel, + Transport* TransportFactory::newTransport(std::string transportModel, thermo_t* phase, int log_level) { if (transportModel == "") return new Transport; @@ -342,6 +336,11 @@ namespace Cantera { dtr = (DustyGasTransport*)tr; dtr->initialize(phase, gastr); break; + case cSimpleTransport: + tr = new SimpleTransport(); + initLiquidTransport(tr, phase, log_level); + tr->setThermo(*phase); + break; #ifdef WITH_IDEAL_SOLUTIONS case cLiquidTransport: tr = new LiquidTransport; @@ -357,27 +356,48 @@ namespace Cantera { break; #endif default: - throw CanteraError("newTransport","unknown transport model"); + throw CanteraError("newTransport","unknown transport model: " + transportModel); } phase->restoreState(state); return tr; } + /** + * make one of several transport models, and return a base class + * pointer to it. + */ + Transport* TransportFactory::newTransport(thermo_t* phase, int log_level) { + XML_Node &phaseNode=phase->xml(); + /* + * Find the Thermo XML node + */ + if (!phaseNode.hasChild("transport")) { + throw CanteraError("TransportFactory::newTransport", + "no transport XML node"); + } + XML_Node& transportNode = phaseNode.child("transport"); + string transportModel = transportNode.attrib("model"); + if (transportModel == "") { + throw CanteraError("TransportFactory::newTransport", + "transport XML node doesn't have a model string"); + } + return newTransport(transportModel, phase,log_level); + } -/** + /** * Prepare to build a new kinetic-theory-based transport manager * for low-density gases. Uses polynomial fits to Monchick & Mason * collision integrals. */ void TransportFactory::setupMM(std::ostream &flog, const std::vector &transport_database, - thermo_t* thermo, int mode, int log_level, TransportParams& tr) { + thermo_t* thermo, int mode, int log_level, GasTransportParams& tr) { // constant mixture attributes tr.thermo = thermo; - tr.nsp = tr.thermo->nSpecies(); - int nsp = tr.nsp; + tr.nsp_ = tr.thermo->nSpecies(); + int nsp = tr.nsp_; tr.tmin = thermo->minTemp(); tr.tmax = thermo->maxTemp(); @@ -387,7 +407,7 @@ namespace Cantera { copy(tr.thermo->molecularWeights().begin(), tr.thermo->molecularWeights().end(), tr.mw.begin()); - tr.mode = mode; + tr.mode_ = mode; tr.epsilon.resize(nsp, nsp, 0.0); tr.delta.resize(nsp, nsp, 0.0); tr.reducedMass.resize(nsp, nsp, 0.0); @@ -456,7 +476,6 @@ namespace Cantera { } } - // Chemkin fits the entire T* range in the Monchick and Mason tables, // so modify tstar_min and tstar_max if in Chemkin compatibility mode @@ -496,30 +515,155 @@ namespace Cantera { } + + /** + * Prepare to build a new transport manager for liquids assuming that + * viscosity transport data is provided in Arhennius form. + */ + void TransportFactory::setupLiquidTransport(std::ostream &flog, + const std::vector &transport_database, + thermo_t* thermo, int log_level, LiquidTransportParams& trParam) { + + // constant mixture attributes + trParam.thermo = thermo; + trParam.nsp_ = trParam.thermo->nSpecies(); + int nsp = trParam.nsp_; + + trParam.tmin = thermo->minTemp(); + trParam.tmax = thermo->maxTemp(); + trParam.mw.resize(nsp); + trParam.log_level = log_level; + + copy(trParam.thermo->molecularWeights().begin(), + trParam.thermo->molecularWeights().end(), trParam.mw.begin()); + + //trParam.epsilon.resize(nsp, nsp, 0.0); + //trParam.delta.resize(nsp, nsp, 0.0); + //trParam.reducedMass.resize(nsp, nsp, 0.0); + //trParam.dipole.resize(nsp, nsp, 0.0); + //trParam.diam.resize(nsp, nsp, 0.0); + //trParam.polar.resize(nsp, false); + //trParam.poly.resize(nsp); + //trParam.sigma.resize(nsp); + //trParam.eps.resize(nsp); + + XML_Node root, log; + getLiquidTransportData(transport_database, log, + trParam.thermo->speciesNames(), trParam); + + //int i, j; + //for (i = 0; i < nsp; i++) trParam.poly[i].resize(nsp); + + //doublereal ts1, ts2, tstar_min = 1.e8, tstar_max = 0.0; + //doublereal f_eps, f_sigma; + + //DenseMatrix& diam = trParam.diam; + //DenseMatrix& epsilon = trParam.epsilon; + + //for (i = 0; i < nsp; i++) + // { + // for (j = i; j < nsp; j++) + // { + // // the reduced mass + // trParam.reducedMass(i,j) = + // trParam.mw[i] * trParam.mw[j] / (Avogadro * (trParam.mw[i] + trParam.mw[j])); + // + // // hard-sphere diameter for (i,j) collisions + // diam(i,j) = 0.5*(trParam.sigma[i] + trParam.sigma[j]); + // + // // the effective well depth for (i,j) collisions + // epsilon(i,j) = sqrt(trParam.eps[i]*trParam.eps[j]); + // + // // The polynomial fits of collision integrals vs. T* + // // will be done for the T* from tstar_min to tstar_max + // ts1 = Boltzmann * trParam.tmin/epsilon(i,j); + // ts2 = Boltzmann * trParam.tmax/epsilon(i,j); + // if (ts1 < tstar_min) tstar_min = ts1; + // if (ts2 > tstar_max) tstar_max = ts2; + // + // // the effective dipole moment for (i,j) collisions + // trParam.dipole(i,j) = sqrt(trParam.dipole(i,i)*trParam.dipole(j,j)); + // + // // reduced dipole moment delta* (nondimensional) + // doublereal d = diam(i,j); + // trParam.delta(i,j) = 0.5 * trParam.dipole(i,j)*trParam.dipole(i,j) + // / (epsilon(i,j) * d * d * d); + // + // makePolarCorrections(i, j, trParam, f_eps, f_sigma); + // trParam.diam(i,j) *= f_sigma; + // epsilon(i,j) *= f_eps; + // + // // properties are symmetric + // trParam.reducedMass(j,i) = trParam.reducedMass(i,j); + // diam(j,i) = diam(i,j); + // epsilon(j,i) = epsilon(i,j); + // trParam.dipole(j,i) = trParam.dipole(i,j); + // trParam.delta(j,i) = trParam.delta(i,j); + // } + // } + + // Chemkin fits the entire T* range in the Monchick and Mason tables, + // so modify tstar_min and tstar_max if in Chemkin compatibility mode + + //if (mode == CK_Mode) { + // tstar_min = 0.101; + // tstar_max = 99.9; + //} + + + // initialize the collision integral calculator for the desired + // T* range + //#ifdef DEBUG_MODE + // if (m_verbose) { + // trParam.xml->XML_open(flog, "collision_integrals"); + // } + //#endif + // m_integrals = new MMCollisionInt; + // m_integrals->init(trParam.xml, tstar_min, tstar_max, log_level); + // fitCollisionIntegrals(flog, trParam); + //#ifdef DEBUG_MODE + // if (m_verbose) { + // trParam.xml->XML_close(flog, "collision_integrals"); + // } + //#endif + // // make polynomial fits + //#ifdef DEBUG_MODE + // if (m_verbose) { + // trParam.xml->XML_open(flog, "property fits"); + // } + //#endif + // fitProperties(trParam, flog); + //#ifdef DEBUG_MODE + // if (m_verbose) { + // trParam.xml->XML_close(flog, "property fits"); + // } + //#endif + } + + void TransportFactory::initTransport(Transport* tran, thermo_t* thermo, int mode, int log_level) { const std::vector & transport_database = thermo->speciesData(); - TransportParams tr; + GasTransportParams trParam; #ifdef DEBUG_MODE ofstream flog("transport_log.xml"); - tr.xml = new XML_Writer(flog); + trParam.xml = new XML_Writer(flog); if (m_verbose) { - tr.xml->XML_open(flog, "transport"); + trParam.xml->XML_open(flog, "transport"); } #else // create the object, but don't associate it with a file std::ostream &flog(std::cout); #endif // set up Monchick and Mason collision integrals - setupMM(flog, transport_database, thermo, mode, log_level, tr); - + setupMM(flog, transport_database, thermo, mode, log_level, trParam); // do model-specific initialization - tran->init(tr); + tran->initGas(trParam); #ifdef DEBUG_MODE if (m_verbose) { - tr.xml->XML_close(flog, "transport"); + trParam.xml->XML_close(flog, "transport"); } // finished with log file flog.close(); @@ -528,11 +672,37 @@ namespace Cantera { } - void - TransportFactory::initLiquidTransport(Transport* tran, - thermo_t* thermo, - int log_level) { + /** Similar to initTransport except uses LiquidTransportParams + * class and calls setupLiquidTransport(). + */ + void TransportFactory::initLiquidTransport(Transport* tran, + thermo_t* thermo, + int log_level) { + const std::vector & transport_database = thermo->speciesData(); + + LiquidTransportParams trParam; +#ifdef DEBUG_MODE + ofstream flog("transport_log.xml"); + trParam.xml = new XML_Writer(flog); + if (m_verbose) { + trParam.xml->XML_open(flog, "transport"); + } +#else + // create the object, but don't associate it with a file + std::ostream &flog(std::cout); +#endif + setupLiquidTransport(flog, transport_database, thermo, log_level, trParam); + // do model-specific initialization + tran->initLiquid(trParam); +#ifdef DEBUG_MODE + if (m_verbose) { + trParam.xml->XML_close(flog, "transport"); + } + // finished with log file + flog.close(); +#endif + return; } @@ -546,12 +716,12 @@ namespace Cantera { void TransportFactory::fitCollisionIntegrals(ostream& logfile, - TransportParams& tr) { + GasTransportParams& tr) { vector_fp::iterator dptr; doublereal dstar; - int nsp = tr.nsp; - int mode = tr.mode; + int nsp = tr.nsp_; + int mode = tr.mode_; int i, j; // Chemkin fits to sixth order polynomials @@ -630,7 +800,7 @@ namespace Cantera { * these species read from the file. */ void TransportFactory::getTransportData(const std::vector &xspecies, - XML_Node& log, const std::vector &names, TransportParams& tr) + XML_Node& log, const std::vector &names, GasTransportParams& tr) { string name; int geom; @@ -653,6 +823,7 @@ namespace Cantera { for (i = 0; i < nsp; i++) { const XML_Node& sp = *xspecies[i]; name = sp["name"]; + // std::cout << "Processing node for " << name << std::endl; // put in a try block so that species with no 'transport' // child are skipped, instead of throwing an exception. @@ -698,7 +869,7 @@ namespace Cantera { } } - for (i = 0; i < tr.nsp; i++) { + for (i = 0; i < tr.nsp_; i++) { GasTransportData& trdat = datatable[names[i]]; @@ -743,6 +914,267 @@ namespace Cantera { } } + /** + * Read transport property data from a file for a list of species. + * Given the name of a file containing transport property + * parameters and a list of species names, this method returns an + * instance of TransportParams containing the transport data for + * these species read from the file. + */ + void TransportFactory::getLiquidTransportData( const std::vector &xspecies, + XML_Node& log, + const std::vector &names, + LiquidTransportParams& trParam) + { + std::string name; + /* + * Create a map of species names versus liquid transport data parameters + */ + std::map datatable; + doublereal A_visc, n_visc, Tact_visc, hydrodynamic_radius; + doublereal A_thcond, n_thcond, Tact_thcond; + doublereal A_spdiff, n_spdiff, Tact_spdiff; + + int nsp = static_cast(xspecies.size()); + std::cout << "Size of xspecies " << nsp << std::endl; + + // read all entries in database into 'datatable' and check for + // errors. Note that this procedure validates all entries, not + // only those for the species listed in 'names'. + + int linenum = 0; + int i; + for (i = 0; i < nsp; i++) { + const XML_Node& sp = *xspecies[i]; + name = sp["name"]; + vector_fp vCoeff; + // std::cout << "Processing node for " << name << std::endl; + + // put in a try block so that species with no 'transport' + // child are skipped, instead of throwing an exception. + try { + if (sp.hasChild("transport")) { + XML_Node& trNode = sp.child("transport"); + + // Fill datatable with LiquidTransportData objects for error checking + // and then insertion into LiquidTransportData objects below. + LiquidTransportData data; + data.speciesName = name; + + /* + * hydrodynamic radius + * + * format: + * 3.0 + * 3.0 + */ + if (trNode.hasChild("hydrodynamic_radius")) { + XML_Node& hnode = trNode.child("hydrodynamic_radius"); + std::string model = lowercase(hnode["model"]); + if (model == "" || model == "constant") { + hydrodynamic_radius = hnode.fp_value(); + if (hydrodynamic_radius > 0.0) data.hydroradius = hydrodynamic_radius; + else throw TransportDBError(linenum, + "negative or zero hydrodynamic radius"); + data.model_hydroradius = LTR_MODEL_CONSTANT; + } else { + throw CanteraError(" TransportFactory::getLiquidTransportData", + "Unknown model for hydrodynamic_radius:" + model); + } + } + + /* + * viscosity + * + * format: + * 3.0 + * 3.0 + * + * 1.0 + * 2.0 + * 3.0 + * + * + * + * 0.0. 1.0, 2.0, 3.0, 4.0 + * + * + */ + if (trNode.hasChild("viscosity")) { + XML_Node& vnode = trNode.child("viscosity"); + std::string model = lowercase(vnode["model"]); + if (model == "" || model == "constant") { + A_visc = vnode.fp_value(); + if (A_visc > 0.0) (data.viscCoeffs).push_back(A_visc); + else throw TransportDBError(linenum, + "negative or zero viscosity"); + data.model_viscosity = LTR_MODEL_CONSTANT; + } else if (model == "arrhenius") { + getArrhenius(vnode, A_visc, n_visc, Tact_visc); + if (A_visc <= 0.0) { + throw TransportDBError(linenum, "negative or zero viscosity"); + } + (data.viscCoeffs).push_back(A_visc); + (data.viscCoeffs).push_back(n_visc); + (data.viscCoeffs).push_back(Tact_visc); + data.model_viscosity = LTR_MODEL_ARRHENIUS; + } else if (model == "coeff") { + getFloatArray(vnode, vCoeff, true); + data.viscCoeffs = vCoeff; + vCoeff.clear(); + data.model_viscosity = LTR_MODEL_COEFF; + } else { + throw CanteraError(" TransportFactory::getLiquidTransportData", + "Unknown model for viscosity:" + vnode["model"]); + } + } + + /* + * thermal_conductivity + * + * format: + * 3.0 + * 3.0 + * + * 1.0 + * 2.0 + * 3.0 + * + * + * + * 0.0. 1.0, 2.0, 3.0, 4.0 + * + * + */ + if (trNode.hasChild("thermal_conductivity")) { + XML_Node& tnode = trNode.child("thermal_conductivity"); + std::string model = lowercase(tnode["model"]); + if (model == "" || model == "constant") { + A_thcond = tnode.fp_value(); + if (A_thcond > 0.0) (data.thermalCondCoeffs).push_back(A_thcond); + else throw TransportDBError(linenum, + "negative or zero thermal_conductivity"); + data.model_thermalCond = LTR_MODEL_CONSTANT; + } else if (model == "arrhenius") { + getArrhenius(tnode, A_thcond, n_thcond, Tact_thcond); + if (A_thcond <= 0.0) { + throw TransportDBError(linenum, "negative or zero thermal_conductivity"); + } + (data.thermalCondCoeffs).push_back(A_thcond); + (data.thermalCondCoeffs).push_back(n_thcond); + (data.thermalCondCoeffs).push_back(Tact_thcond); + data.model_thermalCond = LTR_MODEL_ARRHENIUS; + } else if (model == "coeff") { + getFloatArray(tnode, vCoeff, true); + data.thermalCondCoeffs = vCoeff; + vCoeff.clear(); + data.model_thermalCond = LTR_MODEL_COEFF; + } else { + throw CanteraError(" TransportFactory::getLiquidTransportData", + "Unknown model for thermal_conductivity:" + tnode["model"]); + } + } + + + /* + * speciesDiffusivity + * + * format: + * 3.0 + * 3.0 + * + * 1.0 + * 2.0 + * 3.0 + * + * + * + * 0.0. 1.0, 2.0, 3.0, 4.0 + * + * + */ + if (trNode.hasChild("speciesDiffusivity")) { + XML_Node& dnode = trNode.child("speciesDiffusivity"); + std::string model = lowercase(dnode["model"]); + if (model == "" || model == "constant") { + A_spdiff = dnode.fp_value(); + if (A_spdiff > 0.0) (data.speciesDiffusivityCoeffs).push_back(A_spdiff); + else throw TransportDBError(linenum, + "negative or zero speciesDiffusivity"); + data.model_speciesDiffusivity = LTR_MODEL_CONSTANT; + } else if (model == "arrhenius") { + getArrhenius(dnode, A_spdiff, n_spdiff, Tact_spdiff); + if (A_spdiff <= 0.0) { + throw TransportDBError(linenum, "negative or zero speciesDiffusivity"); + } + (data.speciesDiffusivityCoeffs).push_back(A_spdiff); + (data.speciesDiffusivityCoeffs).push_back(n_spdiff); + (data.speciesDiffusivityCoeffs).push_back(Tact_spdiff); + data.model_speciesDiffusivity = LTR_MODEL_ARRHENIUS; + } else if (model == "coeff") { + getFloatArray(dnode, vCoeff, true); + data.speciesDiffusivityCoeffs = vCoeff; + data.model_speciesDiffusivity = LTR_MODEL_COEFF; + } else { + throw CanteraError(" TransportFactory::getLiquidTransportData", + "Unknown model for speciesDiffusivity:" + dnode["model"]); + } + } + + datatable[name] = data; + } + } + catch(CanteraError) { + ; + } + } + + trParam.LTData.clear(); + for (i = 0; i < trParam.nsp_; i++) { + + LiquidTransportData& trdat = datatable[names[i]]; + + // 'datatable' returns a default TransportData object if + // the species name is not one in the transport database. + // This can be detected by examining 'geometry'. + if (trdat.viscCoeffs[0] < 0) { + throw TransportDBError(0,"no transport data found for species " + + names[i]); + } + + // parameters should be converted to SI units before storing + if (trdat.viscCoeffs.size() > 0) { + trParam.visc_A[i] = trdat.viscCoeffs[0] ; + } + if (trdat.viscCoeffs.size() > 2) { + trParam.visc_n[i] = trdat.viscCoeffs[1] ; + trParam.visc_Tact[i] = trdat.viscCoeffs[2] ; + } + + if (trdat.thermalCondCoeffs.size() > 0) { + trParam.thermCond_A[i] = trdat.thermalCondCoeffs[0] ; + } + if (trdat.thermalCondCoeffs.size() > 2) { + trParam.thermCond_n[i] = trdat.thermalCondCoeffs[1] ; + trParam.thermCond_Tact[i] = trdat.thermalCondCoeffs[2] ; + } + + // Angstroms -> meters + trParam.hydroRadius[i] = 1.e-10 * trdat.hydroradius; + + /* + * this is a much more general way to handle the transfer + * -> calling the default copy constructor for LiquidTransportData + */ + trParam.LTData.push_back(trdat); + } + + // Need to identify a method to obtain interaction matrices. + // This will fill LiquidTransportParams members visc_Eij, visc_Sij + trParam.visc_Eij.resize(trParam.nsp_,trParam.nsp_); + cout << "No support for species viscosity interactions in TransportFactory.cpp" << endl; + } + /********************************************************* * @@ -772,7 +1204,7 @@ namespace Cantera { * D(i,j)/sqrt(k_BT)) = \sum_{n = 0}^4 a_n(i,j) (\log T)^n * \f] */ - void TransportFactory::fitProperties(TransportParams& tr, + void TransportFactory::fitProperties(GasTransportParams& tr, ostream& logfile) { doublereal tstar; int k, j, n, ndeg = 0; @@ -782,7 +1214,7 @@ namespace Cantera { // number of points to use in generating fit data const int np = 50; - int mode = tr.mode; + int mode = tr.mode_; int degree = (mode == CK_Mode ? 3 : 4); doublereal t, om22; @@ -834,7 +1266,7 @@ namespace Cantera { c1, cv_rot, cv_int, f_rot, f_trans, om11; doublereal diffcoeff; - for (k = 0; k < tr.nsp; k++) + for (k = 0; k < tr.nsp_; k++) { for (n = 0; n < np; n++) { t = tr.tmin + dt*n; @@ -971,7 +1403,7 @@ namespace Cantera { tr.xml->XML_comment(logfile,s); } if (tr.log_level >= 2) - for (k = 0; k < tr.nsp; k++) { + for (k = 0; k < tr.nsp_; k++) { tr.xml->XML_writeVector(logfile, " ", tr.thermo->speciesName(k), degree+1, DATA_PTR(tr.condcoeffs[k])); } @@ -999,9 +1431,9 @@ namespace Cantera { mxerr = 0.0, mxrelerr = 0.0; vector_fp diff(np + 1); doublereal eps, sigma; - for (k = 0; k < tr.nsp; k++) + for (k = 0; k < tr.nsp_; k++) { - for (j = k; j < tr.nsp; j++) { + for (j = k; j < tr.nsp_; j++) { ipoly = tr.poly[k][j]; for (n = 0; n < np; n++) { diff --git a/Cantera/src/transport/TransportFactory.h b/Cantera/src/transport/TransportFactory.h old mode 100755 new mode 100644 index 879218642..d38603557 --- a/Cantera/src/transport/TransportFactory.h +++ b/Cantera/src/transport/TransportFactory.h @@ -34,6 +34,7 @@ #include "ct_defs.h" #include "TransportBase.h" #include "FactoryBase.h" +#include "LiquidTransportData.h" #if defined(THREAD_SAFE_CANTERA) #include @@ -61,9 +62,11 @@ namespace Cantera { doublereal rotRelaxNumber; }; + // forward references class MMCollisionInt; - class TransportParams; + class GasTransportParams; + class LiquidTransportParams; class XML_Node; @@ -121,23 +124,35 @@ namespace Cantera { * single instance. */ virtual ~TransportFactory(); - - /// Build a new transport manager + //! Build a new transport manager using a transport manager + //! that may not be the same as in the phase description + /*! + * @param model String name for the transport manager + * @param thermo ThermoPhase object + * @param log_level log level + */ virtual Transport* - newTransport(std::string model="", thermo_t* thermo=0, int log_level=0); + newTransport(std::string model, thermo_t* thermo, int log_level=0); + + //! Build a new transport manager using the default transport manager + //! in the phase description + /*! + * @param thermo ThermoPhase object + * @param log_level log level + */ + virtual Transport* + newTransport(thermo_t* thermo, int log_level=0); /// Initialize an existing transport manager virtual void initTransport(Transport* tr, - thermo_t* thermo=0, int mode=0, int log_level=0); + thermo_t* thermo, int mode=0, int log_level=0); /// Initialize an existing transport manager for liquid phase virtual void initLiquidTransport(Transport* tr, - thermo_t* thermo=0, + thermo_t* thermo, int log_level=0); - - private: //! Static instance of the factor -> This is the only instance of this @@ -147,8 +162,6 @@ namespace Cantera { static boost::mutex transport_mutex ; #endif - - //! The constructor is private; use static method factory() to //! get a pointer to a factory instance /*! @@ -161,34 +174,54 @@ namespace Cantera { void getTransportData(const std::vector &db, XML_Node& log, const std::vector& names, - TransportParams& tr); + GasTransportParams& tr); + + + //! Read transport property data from a file for a list of species. + /*! + * + * Given the name of a file containing transport property + * parameters and a list of species names, this method returns an + * instance of TransportParams containing the transport data for + * these species read from the file. + * + */ + void getLiquidTransportData(const std::vector &db, + XML_Node& log, const std::vector& names, + LiquidTransportParams& tr); /** Generate polynomial fits to viscosity, conductivity, and * binary diffusion coefficients */ - void fitProperties(TransportParams& tr, std::ostream & logfile); + void fitProperties(GasTransportParams& tr, std::ostream & logfile); /// Generate polynomial fits to collision integrals void fitCollisionIntegrals(std::ostream & logfile, - TransportParams& tr); + GasTransportParams& tr); void setupMM(std::ostream &flog, const std::vector &transport_database, thermo_t* thermo, int mode, int log_level, - TransportParams& tr); + GasTransportParams& tr); + + + void setupLiquidTransport(std::ostream &flog, const std::vector &transport_database, + thermo_t* thermo, int log_level, + LiquidTransportParams& tr); /// Second-order correction to the binary diffusion coefficients void getBinDiffCorrection(doublereal t, - const TransportParams& tr, int k, int j, + const GasTransportParams& tr, int k, int j, doublereal xk, doublereal xj, doublereal& fkj, doublereal& fjk); /// Corrections for polar-nonpolar binary diffusion coefficients void makePolarCorrections(int i, int j, - const TransportParams& tr, doublereal& f_eps, + const GasTransportParams& tr, doublereal& f_eps, doublereal& f_sigma); + //! Boolean indicating whether to turn on verbose printing bool m_verbose; @@ -205,8 +238,8 @@ namespace Cantera { * Create a new transport manager instance. * @ingroup transportProps */ - inline Transport* newTransportMgr(std::string transportModel="", - thermo_t* thermo=0, int loglevel=0, + inline Transport* newTransportMgr(std::string transportModel = "", + thermo_t* thermo = 0, int loglevel=0, TransportFactory* f=0) { if (f == 0) { f = TransportFactory::factory(); @@ -218,7 +251,25 @@ namespace Cantera { * the need for multiple cantera and transport library statements * for applications that don't have transport in them. */ - //TransportFactory::deleteFactory(); + return ptr; + } + + /** + * Create a new transport manager instance. + * @ingroup transportProps + */ + inline Transport* newDefaultTransportMgr(thermo_t* thermo, int loglevel=0, + TransportFactory* f=0) { + if (f == 0) { + f = TransportFactory::factory(); + } + Transport* ptr = f->newTransport(thermo, loglevel); + /* + * Note: We delete the static s_factory instance here, instead of in + * appdelete() in misc.cpp, to avoid linking problems involving + * the need for multiple cantera and transport library statements + * for applications that don't have transport in them. + */ return ptr; } diff --git a/Cantera/src/transport/TransportParams.h b/Cantera/src/transport/TransportParams.h index eb46624a0..24eeaa1cd 100755 --- a/Cantera/src/transport/TransportParams.h +++ b/Cantera/src/transport/TransportParams.h @@ -3,7 +3,6 @@ #include - #include "ct_defs.h" #include "TransportBase.h" #include "xml.h" @@ -11,10 +10,9 @@ namespace Cantera { - /** - * - * Holds transport data. Used by TransportFactory. - * + /** + * Base class to hold transport model parameters. + * Used by TransportFactory. */ class TransportParams { @@ -22,18 +20,43 @@ namespace Cantera { TransportParams() : thermo(0), xml(0) {} virtual ~TransportParams(); - int nsp; + int nsp_; // phase_t* mix; thermo_t* thermo; vector_fp mw; // polynomial fits - std::vector visccoeffs; - std::vector condcoeffs; - std::vector diffcoeffs; - vector_fp polytempvec; - + //temperature-fit viscosity + std::vector visccoeffs; + //temperature-fit heat conduction + std::vector condcoeffs; + //temperature-fit diffusivity + std::vector diffcoeffs; + vector_fp polytempvec; + + //minimum and maximum temperatures for parameter fits + doublereal tmax, tmin; + int mode_; + XML_Writer* xml; + int log_level; + + }; + + + /** + * Holds transport model parameters relevant to transport in ideal + * gases with a kinetic theory of gases derived transport model. + * Used by TransportFactory. + */ + class GasTransportParams : public TransportParams { + + public: + + GasTransportParams() {} + ~GasTransportParams() {} + + std::vector > poly; std::vector omega22_poly; std::vector astar_poly; @@ -53,11 +76,9 @@ namespace Cantera { DenseMatrix epsilon; DenseMatrix dipole; DenseMatrix delta; - doublereal tmax, tmin; - int mode; - XML_Writer* xml; - int log_level; + }; + } -#endif +#endif //CT_TRANSPORTPARAMS_H diff --git a/Cantera/src/transport/WaterTransport.h b/Cantera/src/transport/WaterTransport.h index 6e3b10b4a..31d89966a 100644 --- a/Cantera/src/transport/WaterTransport.h +++ b/Cantera/src/transport/WaterTransport.h @@ -81,7 +81,7 @@ namespace Cantera { virtual ~WaterTransport(); //! Return the model id for this transport parameterization - virtual int model() { + virtual int model() const { return cWaterTransport; } diff --git a/configure b/configure index b8e9c815c..9473d7c46 100755 --- a/configure +++ b/configure @@ -10672,7 +10672,7 @@ fi - ac_config_files="$ac_config_files Makefile Cantera/Makefile Cantera/src/Makefile Cantera/src/base/Makefile Cantera/src/zeroD/Makefile Cantera/src/oneD/Makefile Cantera/src/converters/Makefile Cantera/src/transport/Makefile Cantera/src/thermo/Makefile Cantera/src/kinetics/Makefile Cantera/src/numerics/Makefile Cantera/src/spectra/Makefile Cantera/src/equil/Makefile Cantera/clib/src/Makefile Cantera/fortran/src/Makefile Cantera/fortran/f77demos/f77demos.mak Cantera/fortran/f77demos/Makefile Cantera/matlab/Makefile Cantera/matlab/setup_matlab.py Cantera/python/Makefile Cantera/python/setup.py Cantera/cxx/Makefile Cantera/cxx/src/Makefile Cantera/cxx/demos/Makefile Cantera/cxx/demos/combustor/Makefile Cantera/cxx/demos/combustor/Makefile.install Cantera/cxx/demos/flamespeed/Makefile Cantera/cxx/demos/flamespeed/Makefile.install Cantera/cxx/demos/kinetics1/Makefile Cantera/cxx/demos/kinetics1/Makefile.install Cantera/cxx/demos/NASA_coeffs/Makefile Cantera/cxx/demos/NASA_coeffs/Makefile.install Cantera/cxx/demos/rankine/Makefile Cantera/cxx/demos/rankine/Makefile.install Cantera/cxx/include/Cantera.mak Cantera/cxx/include/Cantera_bt.mak Cantera/user/Makefile Cantera/python/src/Makefile Cantera/python/examples/Makefile Cantera/python/examples/equilibrium/Makefile Cantera/python/examples/equilibrium/adiabatic_flame/Makefile Cantera/python/examples/equilibrium/multiphase_plasma/Makefile Cantera/python/examples/equilibrium/simple_test/Makefile Cantera/python/examples/equilibrium/stoich_flame/Makefile Cantera/python/examples/gasdynamics/isentropic/Makefile Cantera/python/examples/gasdynamics/soundSpeed/Makefile Cantera/python/examples/flames/adiabatic_flame/Makefile Cantera/python/examples/flames/flame1/Makefile Cantera/python/examples/flames/flame2/Makefile Cantera/python/examples/flames/flame_fixed_T/Makefile Cantera/python/examples/flames/free_h2_air/Makefile Cantera/python/examples/flames/npflame1/Makefile Cantera/python/examples/flames/stflame1/Makefile Cantera/python/examples/fuel_cells/Makefile Cantera/python/examples/liquid_vapor/critProperties/Makefile Cantera/python/examples/liquid_vapor/rankine/Makefile Cantera/python/examples/kinetics/Makefile Cantera/python/examples/misc/Makefile Cantera/python/examples/reactors/combustor_sim/Makefile Cantera/python/examples/reactors/functors_sim/Makefile Cantera/python/examples/reactors/mix1_sim/Makefile Cantera/python/examples/reactors/mix2_sim/Makefile Cantera/python/examples/reactors/piston_sim/Makefile Cantera/python/examples/reactors/reactor1_sim/Makefile Cantera/python/examples/reactors/reactor2_sim/Makefile Cantera/python/examples/reactors/sensitivity_sim/Makefile Cantera/python/examples/reactors/surf_pfr_sim/Makefile Cantera/python/examples/surface_chemistry/diamond_cvd/Makefile Cantera/python/examples/surface_chemistry/catcomb_stagflow/Makefile Cantera/python/examples/transport/Makefile Cantera/python/examples/flames/Makefile Cantera/python/examples/gasdynamics/Makefile Cantera/python/examples/liquid_vapor/Makefile Cantera/python/examples/reactors/Makefile Cantera/python/examples/surface_chemistry/Makefile ext/lapack/Makefile ext/blas/Makefile ext/cvode/Makefile ext/math/Makefile ext/recipes/Makefile ext/tpx/Makefile ext/Makefile ext/f2c_libs/Makefile ext/f2c_blas/Makefile ext/f2c_lapack/Makefile ext/f2c_math/Makefile examples/Makefile examples/cxx/Makefile tools/Makefile tools/doc/Cantera.cfg tools/doc/Makefile tools/src/Makefile tools/src/sample.mak tools/src/finish_install.py tools/src/package4mac tools/templates/f77/demo.mak tools/templates/f90/demo.mak tools/templates/cxx/demo.mak tools/testtools/Makefile data/inputs/Makefile data/inputs/mkxml test_problems/Makefile test_problems/cxx_ex/Makefile test_problems/silane_equil/Makefile test_problems/surfkin/Makefile test_problems/spectroscopy/Makefile test_problems/surfSolverTest/Makefile test_problems/diamondSurf/Makefile test_problems/diamondSurf_dupl/Makefile test_problems/ChemEquil_gri_matrix/Makefile test_problems/ChemEquil_gri_pairs/Makefile test_problems/ChemEquil_ionizedGas/Makefile test_problems/ChemEquil_red1/Makefile test_problems/CpJump/Makefile test_problems/mixGasTransport/Makefile test_problems/multiGasTransport/Makefile test_problems/printUtilUnitTest/Makefile test_problems/fracCoeff/Makefile test_problems/negATest/Makefile test_problems/NASA9poly_test/Makefile test_problems/ck2cti_test/Makefile test_problems/ck2cti_test/runtest test_problems/nasa9_reader/Makefile test_problems/nasa9_reader/runtest test_problems/min_python/Makefile test_problems/min_python/minDiamond/Makefile test_problems/min_python/negATest/Makefile test_problems/pureFluidTest/Makefile test_problems/rankine_democxx/Makefile test_problems/python/Makefile test_problems/cathermo/Makefile test_problems/cathermo/issp/Makefile test_problems/cathermo/ims/Makefile test_problems/cathermo/stoichSubSSTP/Makefile test_problems/cathermo/testIAPWS/Makefile test_problems/cathermo/testIAPWSPres/Makefile test_problems/cathermo/testIAPWSTripP/Makefile test_problems/cathermo/testWaterPDSS/Makefile test_problems/cathermo/testWaterTP/Makefile test_problems/cathermo/HMW_test_1/Makefile test_problems/cathermo/HMW_test_3/Makefile test_problems/cathermo/HMW_graph_GvT/Makefile test_problems/cathermo/HMW_graph_GvI/Makefile test_problems/cathermo/HMW_graph_HvT/Makefile test_problems/cathermo/HMW_graph_CpvT/Makefile test_problems/cathermo/HMW_graph_VvT/Makefile test_problems/cathermo/DH_graph_1/Makefile test_problems/cathermo/DH_graph_acommon/Makefile test_problems/cathermo/DH_graph_NM/Makefile test_problems/cathermo/DH_graph_Pitzer/Makefile test_problems/cathermo/DH_graph_bdotak/Makefile test_problems/cathermo/HMW_dupl_test/Makefile test_problems/cathermo/VPissp/Makefile test_problems/VCSnonideal/Makefile test_problems/VPsilane_test/Makefile test_problems/VPsilane_test/runtest test_problems/VCSnonideal/NaCl_equil/Makefile bin/install_tsc" + ac_config_files="$ac_config_files Makefile Cantera/Makefile Cantera/src/Makefile Cantera/src/base/Makefile Cantera/src/zeroD/Makefile Cantera/src/oneD/Makefile Cantera/src/converters/Makefile Cantera/src/transport/Makefile Cantera/src/thermo/Makefile Cantera/src/kinetics/Makefile Cantera/src/numerics/Makefile Cantera/src/spectra/Makefile Cantera/src/equil/Makefile Cantera/clib/src/Makefile Cantera/fortran/src/Makefile Cantera/fortran/f77demos/f77demos.mak Cantera/fortran/f77demos/Makefile Cantera/matlab/Makefile Cantera/matlab/setup_matlab.py Cantera/python/Makefile Cantera/python/setup.py Cantera/cxx/Makefile Cantera/cxx/src/Makefile Cantera/cxx/demos/Makefile Cantera/cxx/demos/combustor/Makefile Cantera/cxx/demos/combustor/Makefile.install Cantera/cxx/demos/flamespeed/Makefile Cantera/cxx/demos/flamespeed/Makefile.install Cantera/cxx/demos/kinetics1/Makefile Cantera/cxx/demos/kinetics1/Makefile.install Cantera/cxx/demos/NASA_coeffs/Makefile Cantera/cxx/demos/NASA_coeffs/Makefile.install Cantera/cxx/demos/rankine/Makefile Cantera/cxx/demos/rankine/Makefile.install Cantera/cxx/include/Cantera.mak Cantera/cxx/include/Cantera_bt.mak Cantera/user/Makefile Cantera/python/src/Makefile Cantera/python/examples/Makefile Cantera/python/examples/equilibrium/Makefile Cantera/python/examples/equilibrium/adiabatic_flame/Makefile Cantera/python/examples/equilibrium/multiphase_plasma/Makefile Cantera/python/examples/equilibrium/simple_test/Makefile Cantera/python/examples/equilibrium/stoich_flame/Makefile Cantera/python/examples/gasdynamics/isentropic/Makefile Cantera/python/examples/gasdynamics/soundSpeed/Makefile Cantera/python/examples/flames/adiabatic_flame/Makefile Cantera/python/examples/flames/flame1/Makefile Cantera/python/examples/flames/flame2/Makefile Cantera/python/examples/flames/flame_fixed_T/Makefile Cantera/python/examples/flames/free_h2_air/Makefile Cantera/python/examples/flames/npflame1/Makefile Cantera/python/examples/flames/stflame1/Makefile Cantera/python/examples/fuel_cells/Makefile Cantera/python/examples/liquid_vapor/critProperties/Makefile Cantera/python/examples/liquid_vapor/rankine/Makefile Cantera/python/examples/kinetics/Makefile Cantera/python/examples/misc/Makefile Cantera/python/examples/reactors/combustor_sim/Makefile Cantera/python/examples/reactors/functors_sim/Makefile Cantera/python/examples/reactors/mix1_sim/Makefile Cantera/python/examples/reactors/mix2_sim/Makefile Cantera/python/examples/reactors/piston_sim/Makefile Cantera/python/examples/reactors/reactor1_sim/Makefile Cantera/python/examples/reactors/reactor2_sim/Makefile Cantera/python/examples/reactors/sensitivity_sim/Makefile Cantera/python/examples/reactors/surf_pfr_sim/Makefile Cantera/python/examples/surface_chemistry/diamond_cvd/Makefile Cantera/python/examples/surface_chemistry/catcomb_stagflow/Makefile Cantera/python/examples/transport/Makefile Cantera/python/examples/flames/Makefile Cantera/python/examples/gasdynamics/Makefile Cantera/python/examples/liquid_vapor/Makefile Cantera/python/examples/reactors/Makefile Cantera/python/examples/surface_chemistry/Makefile ext/lapack/Makefile ext/blas/Makefile ext/cvode/Makefile ext/math/Makefile ext/recipes/Makefile ext/tpx/Makefile ext/Makefile ext/f2c_libs/Makefile ext/f2c_blas/Makefile ext/f2c_lapack/Makefile ext/f2c_math/Makefile examples/Makefile examples/cxx/Makefile tools/Makefile tools/doc/Cantera.cfg tools/doc/Makefile tools/src/Makefile tools/src/sample.mak tools/src/finish_install.py tools/src/package4mac tools/templates/f77/demo.mak tools/templates/f90/demo.mak tools/templates/cxx/demo.mak tools/testtools/Makefile data/inputs/Makefile data/inputs/mkxml test_problems/Makefile test_problems/cxx_ex/Makefile test_problems/silane_equil/Makefile test_problems/surfkin/Makefile test_problems/spectroscopy/Makefile test_problems/surfSolverTest/Makefile test_problems/diamondSurf/Makefile test_problems/diamondSurf_dupl/Makefile test_problems/ChemEquil_gri_matrix/Makefile test_problems/ChemEquil_gri_pairs/Makefile test_problems/ChemEquil_ionizedGas/Makefile test_problems/ChemEquil_red1/Makefile test_problems/CpJump/Makefile test_problems/mixGasTransport/Makefile test_problems/multiGasTransport/Makefile test_problems/printUtilUnitTest/Makefile test_problems/fracCoeff/Makefile test_problems/negATest/Makefile test_problems/NASA9poly_test/Makefile test_problems/ck2cti_test/Makefile test_problems/ck2cti_test/runtest test_problems/nasa9_reader/Makefile test_problems/nasa9_reader/runtest test_problems/min_python/Makefile test_problems/min_python/minDiamond/Makefile test_problems/min_python/negATest/Makefile test_problems/pureFluidTest/Makefile test_problems/rankine_democxx/Makefile test_problems/python/Makefile test_problems/cathermo/Makefile test_problems/cathermo/issp/Makefile test_problems/cathermo/ims/Makefile test_problems/cathermo/stoichSubSSTP/Makefile test_problems/cathermo/testIAPWS/Makefile test_problems/cathermo/testIAPWSPres/Makefile test_problems/cathermo/testIAPWSTripP/Makefile test_problems/cathermo/testWaterPDSS/Makefile test_problems/cathermo/testWaterTP/Makefile test_problems/cathermo/HMW_test_1/Makefile test_problems/cathermo/HMW_test_3/Makefile test_problems/cathermo/HMW_graph_GvT/Makefile test_problems/cathermo/HMW_graph_GvI/Makefile test_problems/cathermo/HMW_graph_HvT/Makefile test_problems/cathermo/HMW_graph_CpvT/Makefile test_problems/cathermo/HMW_graph_VvT/Makefile test_problems/cathermo/DH_graph_1/Makefile test_problems/cathermo/DH_graph_acommon/Makefile test_problems/cathermo/DH_graph_NM/Makefile test_problems/cathermo/DH_graph_Pitzer/Makefile test_problems/cathermo/DH_graph_bdotak/Makefile test_problems/cathermo/HMW_dupl_test/Makefile test_problems/cathermo/VPissp/Makefile test_problems/cathermo/wtWater/Makefile test_problems/VCSnonideal/Makefile test_problems/VPsilane_test/Makefile test_problems/VPsilane_test/runtest test_problems/VCSnonideal/NaCl_equil/Makefile bin/install_tsc" test "x$prefix" = xNONE && prefix=$ac_default_prefix @@ -11298,6 +11298,7 @@ do "test_problems/cathermo/DH_graph_bdotak/Makefile" ) CONFIG_FILES="$CONFIG_FILES test_problems/cathermo/DH_graph_bdotak/Makefile" ;; "test_problems/cathermo/HMW_dupl_test/Makefile" ) CONFIG_FILES="$CONFIG_FILES test_problems/cathermo/HMW_dupl_test/Makefile" ;; "test_problems/cathermo/VPissp/Makefile" ) CONFIG_FILES="$CONFIG_FILES test_problems/cathermo/VPissp/Makefile" ;; + "test_problems/cathermo/wtWater/Makefile" ) CONFIG_FILES="$CONFIG_FILES test_problems/cathermo/wtWater/Makefile" ;; "test_problems/VCSnonideal/Makefile" ) CONFIG_FILES="$CONFIG_FILES test_problems/VCSnonideal/Makefile" ;; "test_problems/VPsilane_test/Makefile" ) CONFIG_FILES="$CONFIG_FILES test_problems/VPsilane_test/Makefile" ;; "test_problems/VPsilane_test/runtest" ) CONFIG_FILES="$CONFIG_FILES test_problems/VPsilane_test/runtest" ;;