From 904dd218f178edcad80f45c7b0b6c53319117287 Mon Sep 17 00:00:00 2001 From: John Hewson Date: Fri, 18 Sep 2009 23:02:59 +0000 Subject: [PATCH 01/24] Creating a change branch of Cantera for development of molten salt phase transport. From 53e1c8fd78517c2660c624ca32e5012f169fd2a3 Mon Sep 17 00:00:00 2001 From: John Hewson Date: Sat, 19 Sep 2009 00:21:46 +0000 Subject: [PATCH 02/24] Added method getSpeciesFluxesES(int ndim, const doublereal* grad_T, int ldx, const doublereal* grad_X, int ldf, const doublereal* grad_Phi, doublereal* fluxes) after passing electrostatic potentialgradient as well as the rest of the gradients. Added initGas( GasTransportParams& tr ) and initLiquid( LiquidTransportParams& tr ) because of the different TranportParams subclasses. There might be a better way to do this. --- Cantera/src/transport/TransportBase.h | 49 +++++++++++++++++++++++++++ 1 file changed, 49 insertions(+) diff --git a/Cantera/src/transport/TransportBase.h b/Cantera/src/transport/TransportBase.h index 0ac33e1a7..777cc8315 100755 --- 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; @@ -221,6 +223,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. @@ -316,6 +351,20 @@ namespace Cantera { virtual bool init(TransportParams& tr) { err("init"); return false; } + /** + * Called by TransportFactory to set parameters. + */ + virtual bool initGas( GasTransportParams& tr ) + { err("init"); return false; } + + /** + * Called by TransportFactory to set parameters. + */ + virtual bool initLiquid( LiquidTransportParams& tr ) + { err("init"); return false; } + + + /** * Set the phase object. From eaeefd4f47d130b4c86ba5e7a12a1c8d88deb28f Mon Sep 17 00:00:00 2001 From: John Hewson Date: Sat, 19 Sep 2009 00:24:47 +0000 Subject: [PATCH 03/24] Added forward declarations references to class LiquidTransportParams; Changed the init(TransportParams& tr) method to initLiquid( LiquidTransportParams& tr ) Removed a number of gas-related transport parameters from the init/initLiquid() method including m_poly, m_diam, m_eps, m_alpha In AqueousTransport.cpp there was a method to return a struct to a GasTransportData. This now returns LiquidTransportData struct, but is not yet working. In LiquidTransport.cpp removed all references to m_polytempvec, a vector of 5 entries of the form - 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; We might want these back soon. This allowed up to comment out the code in the methods that evaluate the temperature dependence of the viscosity, thermal conductivity and diffusivity using polynomials in these variables. --This line, and those below, will be ignored-- M AqueousTransport.h M LiquidTransport.h M AqueousTransport.cpp M LiquidTransport.cpp --- Cantera/src/transport/AqueousTransport.cpp | 21 ++++++--------- Cantera/src/transport/AqueousTransport.h | 10 ++++--- Cantera/src/transport/LiquidTransport.cpp | 31 +++++++++------------- Cantera/src/transport/LiquidTransport.h | 6 ++--- 4 files changed, 29 insertions(+), 39 deletions(-) diff --git a/Cantera/src/transport/AqueousTransport.cpp b/Cantera/src/transport/AqueousTransport.cpp index abb947491..a5aa9e640 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_phi.resize(m_nsp, m_nsp, 0.0); @@ -585,17 +582,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..446701ddb 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 @@ -275,7 +277,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 +288,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/LiquidTransport.cpp b/Cantera/src/transport/LiquidTransport.cpp index 1b13824d0..19219ff35 100644 --- a/Cantera/src/transport/LiquidTransport.cpp +++ b/Cantera/src/transport/LiquidTransport.cpp @@ -174,13 +174,11 @@ namespace Cantera { m_thermo->molecularWeights().end(), m_mw.begin()); // copy polynomials and parameters into local storage - viscCoeffsVector_ = tr.viscCoeffsVector_; + viscCoeffsVector_ = tr.visccoeffs; m_condcoeffs = tr.condcoeffs; //m_diffcoeffs = tr.diffcoeffs; m_mode = tr.mode; - m_diam = tr.diam; - m_eps = tr.eps; m_phi.resize(m_nsp, m_nsp, 0.0); @@ -547,14 +545,6 @@ namespace Cantera { 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; m_diff_temp_ok = false; @@ -690,17 +680,19 @@ 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_cond[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_cond[k] = m_sqrt_t * m_condcoeffs[k]; } } m_cond_temp_ok = true; m_cond_mix_ok = false; + */ } @@ -713,10 +705,11 @@ namespace Cantera { // 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 +718,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 +727,7 @@ namespace Cantera { m_diff_temp_ok = true; m_diff_mix_ok = false; + */ } @@ -756,16 +749,17 @@ namespace Cantera { 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])); + viscSpecies_[k] = exp(viscCoeffsVector_[k]); m_sqvisc[k] = sqrt(viscSpecies_[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]); + m_sqvisc[k] = m_t14 * viscCoeffsVector_[k]; viscSpecies_[k] = (m_sqvisc[k]*m_sqvisc[k]); } } @@ -788,6 +782,7 @@ namespace Cantera { m_visc_temp_ok = true; m_visc_mix_ok = false; + */ } diff --git a/Cantera/src/transport/LiquidTransport.h b/Cantera/src/transport/LiquidTransport.h index 1b131d764..0c7aaa467 100644 --- a/Cantera/src/transport/LiquidTransport.h +++ b/Cantera/src/transport/LiquidTransport.h @@ -25,6 +25,7 @@ using namespace std; // Cantera includes #include "TransportBase.h" #include "DenseMatrix.h" +#include "TransportParams.h" #include "LiquidTransportParams.h" namespace Cantera { @@ -39,7 +40,7 @@ namespace Cantera { - class TransportParams; + class LiquidTransportParams; //! Class LiquidTransport implements mixture-averaged transport @@ -507,9 +508,6 @@ namespace Cantera { */ vector_fp m_cond; - //! Polynomials of the log of the temperature - vector_fp m_polytempvec; - //! State of the mole fraction vector. int m_iStateMF; From 00c36ac25046c7209909836a8f109097be6593c8 Mon Sep 17 00:00:00 2001 From: John Hewson Date: Sat, 19 Sep 2009 00:25:57 +0000 Subject: [PATCH 04/24] Split the TransportParams class into three entities. The base class of the same name holds core components including the thermo object, XML nodes, species names and molecular weights. A derived class GasTransportParams holds kinetic-theory-of-gases-related parameters. LiquidTransportParams in LiquidTransportParams.h now derives from TransportParams class. Removed kinetic-theory-of-gases-related parameters. --- Cantera/src/transport/LiquidTransportParams.h | 51 ++++++++---------- Cantera/src/transport/TransportParams.h | 53 +++++++++++++------ 2 files changed, 58 insertions(+), 46 deletions(-) diff --git a/Cantera/src/transport/LiquidTransportParams.h b/Cantera/src/transport/LiquidTransportParams.h index e55f0bf58..c3021080f 100644 --- a/Cantera/src/transport/LiquidTransportParams.h +++ b/Cantera/src/transport/LiquidTransportParams.h @@ -5,27 +5,38 @@ #include "ct_defs.h" #include "TransportBase.h" +#include "TransportParams.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 + // std::vector viscParams; + vector_fp visc_A; + vector_fp visc_n; + vector_fp visc_Tact; + vector_fp thermCond_A; + vector_fp thermCond_n; + vector_fp thermCond_Tact; + //Hydrodynamic radius of transported molecule + vector_fp hydroRadius; //! Coefficients for the limiting conductivity of ions //! in solution: A_k @@ -47,26 +58,6 @@ namespace Cantera { vector_fp B_k_cond; - // polynomial fits - std::vector viscCoeffsVector_; - std::vector condcoeffs; - std::vector diffcoeffs ; - - - 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/TransportParams.h b/Cantera/src/transport/TransportParams.h index eb46624a0..269be6751 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 From ecfc9a71f5fe5e8de4e6c319ce9a3ea3bcf0eb19 Mon Sep 17 00:00:00 2001 From: John Hewson Date: Sat, 19 Sep 2009 00:26:33 +0000 Subject: [PATCH 05/24] Changed MixTransport::init(TransportParams& tr) to MixTransport::initGas( GasTransportParams& tr ) --- Cantera/src/transport/MixTransport.cpp | 2 +- Cantera/src/transport/MixTransport.h | 4 ++-- Cantera/src/transport/MultiTransport.cpp | 2 +- Cantera/src/transport/MultiTransport.h | 4 ++-- 4 files changed, 6 insertions(+), 6 deletions(-) diff --git a/Cantera/src/transport/MixTransport.cpp b/Cantera/src/transport/MixTransport.cpp index c7fe8ef8a..91f2e89c4 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; diff --git a/Cantera/src/transport/MixTransport.h b/Cantera/src/transport/MixTransport.h index 5fc8d4b96..a7558ca42 100755 --- 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 @@ -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..91d6d3087 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; diff --git a/Cantera/src/transport/MultiTransport.h b/Cantera/src/transport/MultiTransport.h index 78a0f998a..f362ac253 100755 --- a/Cantera/src/transport/MultiTransport.h +++ b/Cantera/src/transport/MultiTransport.h @@ -40,7 +40,7 @@ namespace Cantera { TRANSOLVE_LU }; - class TransportParams; + class GasTransportParams; ///////////////////////////////////////////////////////////// @@ -162,7 +162,7 @@ namespace Cantera { /** * @internal */ - virtual bool init(TransportParams& tr); + virtual bool initGas( GasTransportParams& tr ); /** From a76b99e6ffa1d441d26662227463c9f8b2a85d31 Mon Sep 17 00:00:00 2001 From: John Hewson Date: Sat, 19 Sep 2009 00:35:27 +0000 Subject: [PATCH 06/24] for all of the kinetic-theory-of-gases-related methods supporting the gas-phase transport, we needed to pass the GasTransportParams rather than TransportParams. We made new methods to initialize liquid transport including TransportFactory::setupLiquidTransport() TransportFactory::initLiquidTransport() TransportFactory::getLiquidTransportData() Added new struct LiquidTransportData. Added getArrhenius method copied from kinetics directory to parse Arrhenius form XML. --- Cantera/src/transport/TransportFactory.cpp | 282 +++++++++++++++++++-- Cantera/src/transport/TransportFactory.h | 51 +++- 2 files changed, 309 insertions(+), 24 deletions(-) diff --git a/Cantera/src/transport/TransportFactory.cpp b/Cantera/src/transport/TransportFactory.cpp index fbdfce92f..714565696 100755 --- a/Cantera/src/transport/TransportFactory.cpp +++ b/Cantera/src/transport/TransportFactory.cpp @@ -148,7 +148,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 +220,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 +262,8 @@ 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["User"] = cUserTransport; m_models["None"] = None; //m_models["Radiative"] = cRadiative; @@ -372,7 +374,7 @@ namespace Cantera { */ 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; @@ -456,7 +458,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 +497,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->init(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 +654,37 @@ namespace Cantera { } - void - TransportFactory::initLiquidTransport(Transport* tran, + /** 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->init(trParam); +#ifdef DEBUG_MODE + if (m_verbose) { + trParam.xml->XML_close(flog, "transport"); + } + // finished with log file + flog.close(); +#endif + return; } @@ -546,7 +698,7 @@ namespace Cantera { void TransportFactory::fitCollisionIntegrals(ostream& logfile, - TransportParams& tr) { + GasTransportParams& tr) { vector_fp::iterator dptr; doublereal dstar; @@ -630,7 +782,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; @@ -743,6 +895,104 @@ 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) + { + string name; + std::map datatable; + doublereal A_visc, n_visc, Tact_visc, hydrodynamic_radius; + doublereal A_thcond, n_thcond, Tact_thcond; + + int nsp = static_cast(xspecies.size()); + + // 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"]; + + // put in a try block so that species with no 'transport' + // child are skipped, instead of throwing an exception. + try { + XML_Node& trParam = sp.child("transport"); + + hydrodynamic_radius = getFloat(trParam, "hydrodynamic_radius"); + + XML_Node& visc = trParam.child("viscosity"); + getArrhenius(visc, A_visc, n_visc, Tact_visc ); + + XML_Node& thermCond = trParam.child("thermal_conductivity"); + getArrhenius(thermCond, A_thcond, n_thcond, Tact_thcond ); + + LiquidTransportData data; + data.speciesName = name; + + if ( hydrodynamic_radius > 0.0) data.hydroradius = hydrodynamic_radius; + else throw TransportDBError(linenum, + "negative or zero hydrodynamic radius"); + + if (A_visc >= 0.0) { + data.viscCoeffs[0] = A_visc; + data.viscCoeffs[1] = n_visc; + data.viscCoeffs[2] = Tact_visc; + } + else throw TransportDBError(linenum, + "negative pre-exponential for viscosity"); + + if (A_thcond >= 0.0) { + data.thermalCondCoeffs[0] = A_thcond; + data.thermalCondCoeffs[1] = n_thcond; + data.thermalCondCoeffs[2] = Tact_thcond; + } + else throw TransportDBError(linenum, + "negative pre-exponential for thermalCondoctivity"); + + datatable[name] = data; + } + catch(CanteraError) { + ; + } + } + + 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 + + trParam.visc_A[i] = trdat.viscCoeffs[0] ; + trParam.visc_n[i] = trdat.viscCoeffs[1] ; + trParam.visc_Tact[i] = trdat.viscCoeffs[2] ; + + trParam.thermCond_A[i] = trdat.thermalCondCoeffs[0] ; + trParam.thermCond_n[i] = trdat.thermalCondCoeffs[1] ; + trParam.thermCond_Tact[i] = trdat.thermalCondCoeffs[2] ; + + // Angstroms -> meters + trParam.hydroRadius[i] = 1.e-10 * trdat.hydroradius; + + } + } + /********************************************************* * @@ -772,7 +1022,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; diff --git a/Cantera/src/transport/TransportFactory.h b/Cantera/src/transport/TransportFactory.h index 879218642..489ef7156 100755 --- a/Cantera/src/transport/TransportFactory.h +++ b/Cantera/src/transport/TransportFactory.h @@ -61,10 +61,20 @@ namespace Cantera { doublereal rotRelaxNumber; }; + struct LiquidTransportData { + LiquidTransportData() : speciesName("-"), + hydroradius(-1) {} + std::string speciesName; + doublereal hydroradius; + vector_fp viscCoeffs; + vector_fp thermalCondCoeffs; + }; + // forward references class MMCollisionInt; - class TransportParams; - class XML_Node; + class GasTransportParams; + class LiquidTransportParams; + class XML_Node; //! The purpose of TransportFactory is to create new instances of @@ -161,34 +171,59 @@ namespace Cantera { void getTransportData(const std::vector &db, XML_Node& log, const std::vector& names, - TransportParams& tr); + GasTransportParams& tr); + + 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); + /** + * 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; + } + + //! Boolean indicating whether to turn on verbose printing bool m_verbose; From baa60182c31cc8777282865521f16a983f8d7496 Mon Sep 17 00:00:00 2001 From: John Hewson Date: Tue, 22 Sep 2009 18:22:38 +0000 Subject: [PATCH 07/24] Renaming TransportParams member variables nsp_ and mode_. --- Cantera/src/transport/TransportFactory.cpp | 35 ++++++++++++---------- Cantera/src/transport/TransportParams.h | 2 +- 2 files changed, 20 insertions(+), 17 deletions(-) diff --git a/Cantera/src/transport/TransportFactory.cpp b/Cantera/src/transport/TransportFactory.cpp index 714565696..392f1dcad 100755 --- a/Cantera/src/transport/TransportFactory.cpp +++ b/Cantera/src/transport/TransportFactory.cpp @@ -378,8 +378,8 @@ namespace Cantera { // 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(); @@ -389,7 +389,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); @@ -508,8 +508,8 @@ namespace Cantera { // constant mixture attributes trParam.thermo = thermo; - trParam.nsp = trParam.thermo->nSpecies(); - int nsp = trParam.nsp; + trParam.nsp_ = trParam.thermo->nSpecies(); + int nsp = trParam.nsp_; trParam.tmin = thermo->minTemp(); trParam.tmax = thermo->maxTemp(); @@ -702,8 +702,8 @@ namespace Cantera { 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 @@ -850,7 +850,7 @@ namespace Cantera { } } - for (i = 0; i < tr.nsp; i++) { + for (i = 0; i < tr.nsp_; i++) { GasTransportData& trdat = datatable[names[i]]; @@ -902,8 +902,10 @@ namespace Cantera { * 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) + void TransportFactory::getLiquidTransportData( const std::vector &xspecies, + XML_Node& log, + const std::vector &names, + LiquidTransportParams& trParam) { string name; std::map datatable; @@ -922,6 +924,7 @@ namespace Cantera { const XML_Node& sp = *xspecies[i]; name = sp["name"]; + std::cout << "Processing Liquid Transport for " << name; // put in a try block so that species with no 'transport' // child are skipped, instead of throwing an exception. try { @@ -965,7 +968,7 @@ namespace Cantera { } } - for (i = 0; i < trParam.nsp; i++) { + for (i = 0; i < trParam.nsp_; i++) { LiquidTransportData& trdat = datatable[names[i]]; @@ -1032,7 +1035,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; @@ -1084,7 +1087,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; @@ -1221,7 +1224,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])); } @@ -1249,9 +1252,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/TransportParams.h b/Cantera/src/transport/TransportParams.h index 269be6751..24eeaa1cd 100755 --- a/Cantera/src/transport/TransportParams.h +++ b/Cantera/src/transport/TransportParams.h @@ -37,7 +37,7 @@ namespace Cantera { //minimum and maximum temperatures for parameter fits doublereal tmax, tmin; - int mode; + int mode_; XML_Writer* xml; int log_level; From 330e7666aa20e0c5660d58fc5a02fb8d9dc23df1 Mon Sep 17 00:00:00 2001 From: John Hewson Date: Tue, 22 Sep 2009 18:32:45 +0000 Subject: [PATCH 08/24] Renaming TransportParams member variables nsp_ and mode_. --- Cantera/src/transport/AqueousTransport.cpp | 2 +- Cantera/src/transport/LiquidTransport.cpp | 2 +- Cantera/src/transport/LiquidTransport.h | 5 ++++- Cantera/src/transport/MixTransport.cpp | 2 +- Cantera/src/transport/MultiTransport.cpp | 2 +- 5 files changed, 8 insertions(+), 5 deletions(-) diff --git a/Cantera/src/transport/AqueousTransport.cpp b/Cantera/src/transport/AqueousTransport.cpp index a5aa9e640..171b256a6 100644 --- a/Cantera/src/transport/AqueousTransport.cpp +++ b/Cantera/src/transport/AqueousTransport.cpp @@ -91,7 +91,7 @@ namespace Cantera { m_condcoeffs = tr.condcoeffs; m_diffcoeffs = tr.diffcoeffs; - m_mode = tr.mode; + m_mode = tr.mode_; m_phi.resize(m_nsp, m_nsp, 0.0); diff --git a/Cantera/src/transport/LiquidTransport.cpp b/Cantera/src/transport/LiquidTransport.cpp index 19219ff35..7aede695c 100644 --- a/Cantera/src/transport/LiquidTransport.cpp +++ b/Cantera/src/transport/LiquidTransport.cpp @@ -178,7 +178,7 @@ namespace Cantera { m_condcoeffs = tr.condcoeffs; //m_diffcoeffs = tr.diffcoeffs; - m_mode = tr.mode; + m_mode = tr.mode_; m_phi.resize(m_nsp, m_nsp, 0.0); diff --git a/Cantera/src/transport/LiquidTransport.h b/Cantera/src/transport/LiquidTransport.h index 0c7aaa467..e7aa8dd92 100644 --- a/Cantera/src/transport/LiquidTransport.h +++ b/Cantera/src/transport/LiquidTransport.h @@ -390,7 +390,10 @@ namespace Cantera { * added. */ vector m_diffcoeffs; - + + //! Temperature polynomial for transport property temperature fits. + vector_fp m_polytempvec; + //! Internal value of the gradient of the mole fraction vector /*! diff --git a/Cantera/src/transport/MixTransport.cpp b/Cantera/src/transport/MixTransport.cpp index 91f2e89c4..2dd2da5da 100755 --- a/Cantera/src/transport/MixTransport.cpp +++ b/Cantera/src/transport/MixTransport.cpp @@ -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/MultiTransport.cpp b/Cantera/src/transport/MultiTransport.cpp index 91d6d3087..1f9e443c0 100755 --- a/Cantera/src/transport/MultiTransport.cpp +++ b/Cantera/src/transport/MultiTransport.cpp @@ -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; From 1dad7e320a6abd864d7b75729acd7f3563dbae61 Mon Sep 17 00:00:00 2001 From: John Hewson Date: Thu, 24 Sep 2009 20:46:26 +0000 Subject: [PATCH 09/24] LiquidTransport.h LiquidTransport.cpp Removed a number of variables relevant only to the gas-phase transport coefficient models. These include the following: - vector_fp m_cond; - vector_fp m_polytempvec; - vector m_condcoeffs; - std::vector viscCoeffsVector_; - vector_fp m_sqvisc; - vector_fp viscSpecies_; - DenseMatrix m_wratkj1; - DenseMatrix m_phi; - DenseMatrix m_wratjk; Removed much code related to these variables. Added a number of variables relevant to the liquid-phase transport models including: + vector_fp m_condSpecies; + vector_fp m_viscSpecies; + vector_fp m_visc_A; + vector_fp m_visc_n; + vector_fp m_visc_Tact; + vector_fp m_thermCond_A; + vector_fp m_thermCond_n; + vector_fp m_thermCond_Tact; Changed some of the relevant comments to pertain to the liquid-phase models. --- Cantera/src/transport/LiquidTransport.cpp | 113 ++++++----------- Cantera/src/transport/LiquidTransport.h | 148 +++++++++------------- 2 files changed, 104 insertions(+), 157 deletions(-) diff --git a/Cantera/src/transport/LiquidTransport.cpp b/Cantera/src/transport/LiquidTransport.cpp index 7aede695c..bf69d190e 100644 --- a/Cantera/src/transport/LiquidTransport.cpp +++ b/Cantera/src/transport/LiquidTransport.cpp @@ -100,27 +100,26 @@ 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_n = right.m_visc_n; + m_visc_Tact = right.m_visc_Tact; + m_thermCond_A = right.m_thermCond_A; + m_thermCond_n = right.m_thermCond_n; + m_thermCond_Tact = right.m_thermCond_Tact; 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_viscSpecies = right.m_viscSpecies; m_cond = right.m_cond; - m_polytempvec = right.m_polytempvec; 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; @@ -173,30 +172,29 @@ namespace Cantera { copy(m_thermo->molecularWeights().begin(), m_thermo->molecularWeights().end(), m_mw.begin()); - // copy polynomials and parameters into local storage - viscCoeffsVector_ = tr.visccoeffs; - 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 ; + + m_thermCond_A = tr.thermCond_A ; + m_thermCond_n = tr.thermCond_n ; + m_thermCond_Tact = tr.thermCond_Tact ; + //m_diffcoeffs = tr.diffcoeffs; m_mode = tr.mode_; - m_phi.resize(m_nsp, m_nsp, 0.0); + m_visc_A.resize(m_nsp); + m_visc_n.resize(m_nsp); + m_visc_Tact.resize(m_nsp); + m_thermCond_A.resize(m_nsp); + m_thermCond_n.resize(m_nsp); + m_thermCond_Tact.resize(m_nsp); - 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_condSpecies.resize(m_nsp); m_bdiff.resize(m_nsp, m_nsp); m_molefracs.resize(m_nsp); @@ -247,7 +245,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(); } @@ -256,16 +254,18 @@ 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; + err("constant viscosity not implemented for LiquidTransport."); + //return m_viscmix; } else if (viscosityModel_ == LVISC_MIXTUREAVG) { - m_viscmix = dot_product(viscSpecies_, m_molefracs); + m_viscmix = dot_product(m_viscSpecies, m_molefracs); + } else if (viscosityModel_ == LVISC_INTERACTION) { + m_viscmix = dot_product(m_viscSpecies, m_molefracs); + //now sum over i,j : Gij*Xi*Xj } 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("Wilkes method not implemented for LiquidTransport."); } return m_viscmix; @@ -276,7 +276,7 @@ namespace Cantera { if (!m_visc_temp_ok) { updateViscosity_temp(); } - copy(viscSpecies_.begin(), viscSpecies_.end(), visc); + copy(m_viscSpecies.begin(), m_viscSpecies.end(), visc); } @@ -354,8 +354,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; @@ -683,11 +683,11 @@ namespace Cantera { /* if (m_mode == CK_Mode) { for (k = 0; k < m_nsp; k++) { - m_cond[k] = exp(m_condcoeffs[k]); + m_condSpecies[k] = exp(m_condcoeffs[k]); } } else { for (k = 0; k < m_nsp; k++) { - m_cond[k] = m_sqrt_t * m_condcoeffs[k]; + m_condSpecies[k] = m_sqrt_t * m_condcoeffs[k]; } } m_cond_temp_ok = true; @@ -747,42 +747,13 @@ 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(viscCoeffsVector_[k]); - m_sqvisc[k] = sqrt(viscSpecies_[k]); - } + 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 ); } - else { - for (k = 0; k < m_nsp; k++) { - // the polynomial fit is done for sqrt(visc/sqrt(T)) - m_sqvisc[k] = m_t14 * 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); - } - } - m_visc_temp_ok = true; m_visc_mix_ok = false; - */ } diff --git a/Cantera/src/transport/LiquidTransport.h b/Cantera/src/transport/LiquidTransport.h index e7aa8dd92..424338fbe 100644 --- a/Cantera/src/transport/LiquidTransport.h +++ b/Cantera/src/transport/LiquidTransport.h @@ -137,11 +137,15 @@ namespace Cantera { public: //! default constructor + /*! + * @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); @@ -149,8 +153,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); @@ -170,6 +174,19 @@ 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() { return cLiquidTransport; @@ -179,17 +196,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 @@ -198,25 +212,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 @@ -232,6 +232,20 @@ namespace Cantera { virtual void getMixDiffCoeffs(doublereal* const d); + 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(); + //! Get the Mobilities /*! * @param mobil @@ -332,20 +346,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(); @@ -369,19 +369,16 @@ 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_n; + vector_fp m_visc_Tact; + + //! Pure species thermal conductivities in Arrhenius temperature-dependent form. + vector_fp m_thermCond_A; + vector_fp m_thermCond_n; + vector_fp m_thermCond_Tact; - //! Polynomial coefficients of the conductivities - /*! - * These express the temperature dependendence of the pures - * species conductivities - */ - vector m_condcoeffs; //! Polynomial coefficients of the binary diffusion coefficients /*! @@ -391,10 +388,6 @@ namespace Cantera { */ vector m_diffcoeffs; - //! Temperature polynomial for transport property temperature fits. - vector_fp m_polytempvec; - - //! Internal value of the gradient of the mole fraction vector /*! * Note, this is the only gradient value that can and perhaps @@ -486,19 +479,7 @@ 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; //! Internal value of the species individual thermal conductivities /*! @@ -509,7 +490,7 @@ namespace Cantera { * * controlling update boolean -> m_cond_temp_ok */ - vector_fp m_cond; + vector_fp m_condSpecies; //! State of the mole fraction vector. int m_iStateMF; @@ -593,21 +574,6 @@ 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; @@ -746,6 +712,16 @@ 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. + */ + doublereal err(std::string msg) const; + + }; } #endif From cc680c2b3c2809767dd552522bf55c288c3ec89f Mon Sep 17 00:00:00 2001 From: John Hewson Date: Thu, 24 Sep 2009 20:58:28 +0000 Subject: [PATCH 10/24] Cleanup and debugging of getLiquidTransportData. --- Cantera/src/transport/TransportFactory.cpp | 14 ++++++++------ 1 file changed, 8 insertions(+), 6 deletions(-) diff --git a/Cantera/src/transport/TransportFactory.cpp b/Cantera/src/transport/TransportFactory.cpp index 392f1dcad..19d244fe3 100755 --- a/Cantera/src/transport/TransportFactory.cpp +++ b/Cantera/src/transport/TransportFactory.cpp @@ -805,6 +805,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. @@ -913,7 +914,8 @@ namespace Cantera { doublereal A_thcond, n_thcond, Tact_thcond; 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'. @@ -923,19 +925,19 @@ 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; - std::cout << "Processing Liquid Transport for " << name; // put in a try block so that species with no 'transport' // child are skipped, instead of throwing an exception. try { - XML_Node& trParam = sp.child("transport"); + XML_Node& trNode = sp.child("transport"); - hydrodynamic_radius = getFloat(trParam, "hydrodynamic_radius"); + hydrodynamic_radius = getFloat(trNode, "hydrodynamic_radius"); - XML_Node& visc = trParam.child("viscosity"); + XML_Node& visc = trNode.child("viscosity"); getArrhenius(visc, A_visc, n_visc, Tact_visc ); - XML_Node& thermCond = trParam.child("thermal_conductivity"); + XML_Node& thermCond = trNode.child("thermal_conductivity"); getArrhenius(thermCond, A_thcond, n_thcond, Tact_thcond ); LiquidTransportData data; From 1d1fdcd8f3e2d4af13acb410c5c5c8fad908d6b6 Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Thu, 24 Sep 2009 23:01:26 +0000 Subject: [PATCH 11/24] A little bit of housecleaning. Took out the file Makefile, which shouldn't be in the svn system, and added some ignore files. --- test_problems/cathermo/wtWater/Makefile | 112 ------------------------ 1 file changed, 112 deletions(-) delete mode 100644 test_problems/cathermo/wtWater/Makefile diff --git a/test_problems/cathermo/wtWater/Makefile b/test_problems/cathermo/wtWater/Makefile deleted file mode 100644 index ec802b6a1..000000000 --- a/test_problems/cathermo/wtWater/Makefile +++ /dev/null @@ -1,112 +0,0 @@ -#!/bin/sh - -############################################################################ -# -# Makefile to compile and link a C++ application to -# Cantera. -# -############################################################################# - -# addition to suffixes -.SUFFIXES : .d - -# the name of the executable program to be created -PROG_NAME = wtWater - -# the object files to be linked together. List those generated from Fortran -# and from C/C++ separately -OBJS = wtWater.o - -# additional flags to be passed to the linker. If your program -# requires other external libraries, put them here -LINK_OPTIONS = - -############################################################################# - -# Check to see whether we are in the msvc++ environment -os_is_win = 0 - -# Fortran libraries -ifeq (0, 0) -FORT_LIBS = -else -FORT_LIBS = -endif - -# the C++ compiler -CXX = /home/sntools/extras/compilers/gcc-4.2.4/bin/gcc - -# C++ compile flags -CXX_FLAGS = -g -Wall -DDEBUG_HKM -DDEBUG_HKM_EPEQUIL - -# Ending C++ linking libraries -LCXX_END_LIBS = -lctf2c -lgfortran -lm -lstdc++ - -# the directory where the Cantera libraries are located -CANTERA_LIBDIR=/ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu - -# required Cantera libraries -CANTERA_LIBS = -luser -loneD -lzeroD -lequil -lVCSnonideal -lkinetics -ltransport -lthermo -lctnumerics -lctmath -ltpx -lctspectra -lconverters -lctbase -lsundials_cvodes -lsundials_nvecserial -lctlapack -lctblas -lctf2c -lctcxx - -# Dependencies for CANTERA_LIBS -CANTERA_LIBS_DEP = /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/libuser.a /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/liboneD.a /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/libzeroD.a /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/libequil.a /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/libVCSnonideal.a /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/libkinetics.a /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/libtransport.a /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/libthermo.a /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/libctnumerics.a /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/libctmath.a /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/libtpx.a /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/libctspectra.a /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/libconverters.a /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/libctbase.a /ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu/libctf2c.a $(CANTERA_LIBDIR)/libctcxx.a - -# the directory where Cantera include files may be found. -CANTERA_INCDIR=/ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/include/cantera -INCLUDES=-I$(CANTERA_INCDIR) -I$(CANTERA_INCDIR)/kernel - -# flags passed to the C++ compiler/linker for the linking step -LCXX_FLAGS = -L$(CANTERA_LIBDIR) -L/ascldap/users/hkmoffa/Cantera/gc/canteraDevelop/build/lib/x86_64-unknown-linux-gnu -L/home/hkmoffa/arch/linux64/sundials/lib -g -Wall -DDEBUG_HKM -DDEBUG_HKM_EPEQUIL - -# How to compile C++ source files to object files -.cpp.o: - $(CXX) -c $< $(INCLUDES) $(CXX_FLAGS) - -# How to compile the dependency file -.cpp.d: - /home/sntools/extras/compilers/gcc-4.2.4/bin/g++ -MM $(INCLUDES) $(CXX_FLAGS) $*.cpp > $*.d - -# List of dependency files to be created -DEPENDS=$(OBJS:.o=.d) - -# Program Name -PROGRAM = $(PROG_NAME)$(EXE_EXT) - -all: $(PROGRAM) .depends - -$(PROGRAM): $(OBJS) $(CANTERA_LIBS_DEP) - $(CXX) -o $(PROGRAM) $(OBJS) $(LCXX_FLAGS) $(LINK_OPTIONS) \ - $(CANTERA_LIBS) $(FORT_LIBS) \ - $(LCXX_END_LIBS) - -# Add an additional target for stability: -$(OBJS): $(CANTERA_LIBDIR)/libctbase.a $(CANTERA_LIBDIR)/libthermo.a - -# depends target -> forces recalculation of dependencies -depends: - gmake .depends - -.depends: $(DEPENDS) - cat $(DEPENDS) > .depends - -# Do the test -> For the windows vc++ environment, we have to skip checking on -# whether the program is uptodate, because we don't utilize make -# in that environment to build programs. -test: -ifeq ($(os_is_win), 1) -else - @ gmake -s $(PROGRAM) -endif - @ ./runtest - -clean: - $(RM) $(OBJS) $(PROGRAM) $(DEPENDS) .depends - ../../../bin/rm_cvsignore - (if test -d SunWS_cache ; then \ - $(RM) -rf SunWS_cache ; \ - fi ) - -ifeq ($(wildcard .depends), .depends) -include .depends -endif - From b933e64c790bd9637bbb76dad10ca84dab9bac9d Mon Sep 17 00:00:00 2001 From: John Hewson Date: Fri, 25 Sep 2009 18:40:29 +0000 Subject: [PATCH 12/24] LiquidTransportParams.h Added some comments plus members + DenseMatrix visc_Eij; + DenseMatrix visc_Sij; These members are not populated in any way at this point. A message to this effect is written in TransportFactory::getLiquidTransportData. LiquidTransport.h LiquidTransport.cpp Changed the types of viscosity models from const int LVISC_CONSTANT = 0; - const int LVISC_WILKES = 1; - const int LVISC_MIXTUREAVG = 2; to const int LVISC_CONSTANT = 0; + const int LVISC_INTERACTION = 1; + const int LVISC_AVG_ENERGIES = 2; These are not yet implemented. Added members + DenseMatrix m_visc_Eij; + DenseMatrix m_visc_Sij; + vector_fp m_hydrodynamic_radius; + vector_fp m_visc_logA; //logarithm of coefficient Removed unneeded members: - vector_fp m_eps; - doublereal m_sqrt_t; - doublereal m_t14; - doublereal m_t32; - doublereal m_sqrt_kbt; Changed the name of m_cond to m_condSpecies in agreement with the m_viscSpecies pattern. Implemented some code in LiquidTransport::viscosity() to compute the molecular interactions relevant to mixture viscosity, BUT none of the required parameters or flags are implemented at this point. --- Cantera/src/transport/LiquidTransport.cpp | 77 ++++++++++--------- Cantera/src/transport/LiquidTransport.h | 53 +++++++------ Cantera/src/transport/LiquidTransportParams.h | 27 ++++++- Cantera/src/transport/TransportFactory.cpp | 7 ++ 4 files changed, 103 insertions(+), 61 deletions(-) diff --git a/Cantera/src/transport/LiquidTransport.cpp b/Cantera/src/transport/LiquidTransport.cpp index bf69d190e..aa2f0a3da 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), @@ -101,19 +93,24 @@ namespace Cantera { m_tmax = right.m_tmax; m_mw = right.m_mw; 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; - m_viscSpecies = right.m_viscSpecies; - m_cond = right.m_cond; + 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; @@ -122,13 +119,8 @@ namespace Cantera { viscosityModel_ = right.viscosityModel_; 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; @@ -177,23 +169,28 @@ namespace Cantera { 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_visc_A.resize(m_nsp); - m_visc_n.resize(m_nsp); - m_visc_Tact.resize(m_nsp); - - m_thermCond_A.resize(m_nsp); - m_thermCond_n.resize(m_nsp); - m_thermCond_Tact.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); @@ -257,15 +254,24 @@ namespace Cantera { /* We still need to implement interaction parameters */ /* This constant viscosity model has no input */ if (viscosityModel_ == LVISC_CONSTANT) { + err("constant viscosity not implemented for LiquidTransport."); //return m_viscmix; - } else if (viscosityModel_ == LVISC_MIXTUREAVG) { - m_viscmix = dot_product(m_viscSpecies, m_molefracs); + + } else if (viscosityModel_ == LVISC_AVG_ENERGIES) { + + m_viscmix = exp( dot_product(m_logViscSpecies, m_molefracs) ); + } else if (viscosityModel_ == LVISC_INTERACTION) { - m_viscmix = dot_product(m_viscSpecies, m_molefracs); - //now sum over i,j : Gij*Xi*Xj - } else if (viscosityModel_ == LVISC_WILKES) { - err("Wilkes method not implemented for LiquidTransport."); + + // 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; @@ -540,10 +546,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); // temperature has changed so temp flags are flipped m_visc_temp_ok = false; @@ -749,9 +751,14 @@ namespace Cantera { int k; 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_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] ); } + //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; } diff --git a/Cantera/src/transport/LiquidTransport.h b/Cantera/src/transport/LiquidTransport.h index 424338fbe..06f93bb37 100644 --- a/Cantera/src/transport/LiquidTransport.h +++ b/Cantera/src/transport/LiquidTransport.h @@ -31,8 +31,8 @@ using namespace std; 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; @@ -136,8 +136,13 @@ 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. */ @@ -371,14 +376,32 @@ namespace Cantera { //! 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; + //! 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]. + */ + 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 /*! @@ -469,9 +492,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 @@ -480,6 +503,7 @@ namespace Cantera { * controlling update boolean -> m_visc_temp_ok */ vector_fp m_viscSpecies; + vector_fp m_logViscSpecies; //! Internal value of the species individual thermal conductivities /*! @@ -580,10 +604,6 @@ namespace Cantera { //! 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 @@ -597,21 +617,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; diff --git a/Cantera/src/transport/LiquidTransportParams.h b/Cantera/src/transport/LiquidTransportParams.h index c3021080f..08142b745 100644 --- a/Cantera/src/transport/LiquidTransportParams.h +++ b/Cantera/src/transport/LiquidTransportParams.h @@ -27,14 +27,37 @@ namespace Cantera { //section for liquid transport properties - //Arrhenius parameters for transport coefficients - // std::vector viscParams; + //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; diff --git a/Cantera/src/transport/TransportFactory.cpp b/Cantera/src/transport/TransportFactory.cpp index 19d244fe3..aa7490c2b 100755 --- a/Cantera/src/transport/TransportFactory.cpp +++ b/Cantera/src/transport/TransportFactory.cpp @@ -940,6 +940,8 @@ namespace Cantera { XML_Node& thermCond = trNode.child("thermal_conductivity"); getArrhenius(thermCond, A_thcond, n_thcond, Tact_thcond ); + // Fill datatable with LiquidTransportData objects for error checking + // and then insertion into LiquidTransportData objects below. LiquidTransportData data; data.speciesName = name; @@ -996,6 +998,11 @@ namespace Cantera { trParam.hydroRadius[i] = 1.e-10 * trdat.hydroradius; } + + // 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; } From 26fb3aec3b6f9953c8c96ab95e98a34517525a44 Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Fri, 25 Sep 2009 18:41:20 +0000 Subject: [PATCH 13/24] Updating this directory From 431f167c0f001648ea3ed33fd56bc1448bea70c3 Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Mon, 28 Sep 2009 16:09:58 +0000 Subject: [PATCH 14/24] Changed the method model() to be a const function. Moved Transport methods back into the file TransportBase.cpp Defined err() functions for LiquidTransport object so that the code will link again. --- Cantera/src/transport/AqueousTransport.h | 2 +- Cantera/src/transport/DustyGasTransport.h | 2 +- Cantera/src/transport/LiquidTransport.cpp | 16 ++++++++++ Cantera/src/transport/LiquidTransport.h | 2 +- Cantera/src/transport/MixTransport.h | 2 +- Cantera/src/transport/MultiTransport.h | 2 +- Cantera/src/transport/SolidTransport.h | 2 +- Cantera/src/transport/TransportBase.cpp | 36 ++++++++++++++++++++-- Cantera/src/transport/TransportBase.h | 2 +- Cantera/src/transport/TransportFactory.cpp | 25 --------------- Cantera/src/transport/WaterTransport.h | 2 +- 11 files changed, 58 insertions(+), 35 deletions(-) diff --git a/Cantera/src/transport/AqueousTransport.h b/Cantera/src/transport/AqueousTransport.h index 446701ddb..07e3e3dd9 100644 --- a/Cantera/src/transport/AqueousTransport.h +++ b/Cantera/src/transport/AqueousTransport.h @@ -135,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 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 aa2f0a3da..8175c32eb 100644 --- a/Cantera/src/transport/LiquidTransport.cpp +++ b/Cantera/src/transport/LiquidTransport.cpp @@ -926,4 +926,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 06f93bb37..852641192 100644 --- a/Cantera/src/transport/LiquidTransport.h +++ b/Cantera/src/transport/LiquidTransport.h @@ -193,7 +193,7 @@ namespace Cantera { //! Return the model id for this transport parameterization - virtual int model() { + virtual int model() const { return cLiquidTransport; } diff --git a/Cantera/src/transport/MixTransport.h b/Cantera/src/transport/MixTransport.h index a7558ca42..18a781577 100755 --- a/Cantera/src/transport/MixTransport.h +++ b/Cantera/src/transport/MixTransport.h @@ -51,7 +51,7 @@ namespace Cantera { virtual ~MixTransport() {} - virtual int model() { return cMixtureAveraged; } + virtual int model() const { return cMixtureAveraged; } //! Viscosity of the mixture /*! diff --git a/Cantera/src/transport/MultiTransport.h b/Cantera/src/transport/MultiTransport.h index f362ac253..5979b8ce9 100755 --- a/Cantera/src/transport/MultiTransport.h +++ b/Cantera/src/transport/MultiTransport.h @@ -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 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 index 777cc8315..dd0d7aae3 100755 --- a/Cantera/src/transport/TransportBase.h +++ b/Cantera/src/transport/TransportBase.h @@ -109,7 +109,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 diff --git a/Cantera/src/transport/TransportFactory.cpp b/Cantera/src/transport/TransportFactory.cpp index aa7490c2b..2beffa8e9 100755 --- a/Cantera/src/transport/TransportFactory.cpp +++ b/Cantera/src/transport/TransportFactory.cpp @@ -95,31 +95,6 @@ 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; - } //////////////////// class TransportFactory methods ////////////// 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; } From ff825096006bef0a3aeba164df6e1ad3f08e7de3 Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Fri, 2 Oct 2009 21:08:16 +0000 Subject: [PATCH 15/24] Test commit with small functionality. --- Cantera/src/numerics/Integrator.h | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Cantera/src/numerics/Integrator.h b/Cantera/src/numerics/Integrator.h index 47fc811d8..cd442db88 100755 --- a/Cantera/src/numerics/Integrator.h +++ b/Cantera/src/numerics/Integrator.h @@ -78,8 +78,8 @@ namespace Cantera { virtual void setTolerances(doublereal reltol, int n, doublereal* abstol) { warn("setTolerances"); } - /** - * Set error tolerances. + //! Set error tolerances. + /*! * @param reltol scalar relative tolerance * @param abstol scalar absolute tolerance */ From 825d90987508c622ea1135f01ba543333e1d777c Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Fri, 2 Oct 2009 21:09:12 +0000 Subject: [PATCH 16/24] Took executable bit off of the file. --- Cantera/src/numerics/Integrator.h | 0 1 file changed, 0 insertions(+), 0 deletions(-) mode change 100755 => 100644 Cantera/src/numerics/Integrator.h diff --git a/Cantera/src/numerics/Integrator.h b/Cantera/src/numerics/Integrator.h old mode 100755 new mode 100644 From 636c63482b157fbb73c32083dfe320cb36fb40ad Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Fri, 2 Oct 2009 21:10:26 +0000 Subject: [PATCH 17/24] Took executable bits off --- Cantera/src/kinetics/Enhanced3BConc.h | 0 Cantera/src/kinetics/FalloffFactory.cpp | 0 Cantera/src/kinetics/FalloffFactory.h | 0 Cantera/src/kinetics/FalloffMgr.h | 0 Cantera/src/kinetics/GRI_30_Kinetics.cpp | 0 Cantera/src/kinetics/GRI_30_Kinetics.h | 0 Cantera/src/kinetics/GasKinetics.cpp | 0 Cantera/src/kinetics/GasKinetics.h | 0 Cantera/src/kinetics/GasKineticsWriter.cpp | 0 Cantera/src/kinetics/GasKineticsWriter.h | 0 Cantera/src/kinetics/Group.cpp | 0 Cantera/src/kinetics/Group.h | 0 Cantera/src/kinetics/ImplicitChem.cpp | 0 Cantera/src/kinetics/ImplicitChem.h | 0 Cantera/src/kinetics/ImplicitSurfChem.cpp | 0 Cantera/src/kinetics/ImplicitSurfChem.h | 0 Cantera/src/kinetics/Kinetics.h | 0 Cantera/src/kinetics/RateCoeffMgr.h | 0 Cantera/src/kinetics/ReactionData.h | 0 Cantera/src/kinetics/ReactionPath.cpp | 0 Cantera/src/kinetics/ReactionPath.h | 0 Cantera/src/kinetics/RxnRates.h | 0 Cantera/src/kinetics/StoichManager.h | 0 Cantera/src/kinetics/ThirdBodyMgr.h | 0 Cantera/src/kinetics/reaction_defs.h | 0 25 files changed, 0 insertions(+), 0 deletions(-) mode change 100755 => 100644 Cantera/src/kinetics/Enhanced3BConc.h mode change 100755 => 100644 Cantera/src/kinetics/FalloffFactory.cpp mode change 100755 => 100644 Cantera/src/kinetics/FalloffFactory.h mode change 100755 => 100644 Cantera/src/kinetics/FalloffMgr.h mode change 100755 => 100644 Cantera/src/kinetics/GRI_30_Kinetics.cpp mode change 100755 => 100644 Cantera/src/kinetics/GRI_30_Kinetics.h mode change 100755 => 100644 Cantera/src/kinetics/GasKinetics.cpp mode change 100755 => 100644 Cantera/src/kinetics/GasKinetics.h mode change 100755 => 100644 Cantera/src/kinetics/GasKineticsWriter.cpp mode change 100755 => 100644 Cantera/src/kinetics/GasKineticsWriter.h mode change 100755 => 100644 Cantera/src/kinetics/Group.cpp mode change 100755 => 100644 Cantera/src/kinetics/Group.h mode change 100755 => 100644 Cantera/src/kinetics/ImplicitChem.cpp mode change 100755 => 100644 Cantera/src/kinetics/ImplicitChem.h mode change 100755 => 100644 Cantera/src/kinetics/ImplicitSurfChem.cpp mode change 100755 => 100644 Cantera/src/kinetics/ImplicitSurfChem.h mode change 100755 => 100644 Cantera/src/kinetics/Kinetics.h mode change 100755 => 100644 Cantera/src/kinetics/RateCoeffMgr.h mode change 100755 => 100644 Cantera/src/kinetics/ReactionData.h mode change 100755 => 100644 Cantera/src/kinetics/ReactionPath.cpp mode change 100755 => 100644 Cantera/src/kinetics/ReactionPath.h mode change 100755 => 100644 Cantera/src/kinetics/RxnRates.h mode change 100755 => 100644 Cantera/src/kinetics/StoichManager.h mode change 100755 => 100644 Cantera/src/kinetics/ThirdBodyMgr.h mode change 100755 => 100644 Cantera/src/kinetics/reaction_defs.h 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 From c1d59a912b4daf7b4820d9e49c40b5edd451fb80 Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Sat, 10 Oct 2009 04:23:04 +0000 Subject: [PATCH 18/24] Got most of the way towards creating a new transport Object SimpleTransport that obeys the simplest formulas possible. --- Cantera/src/transport/LiquidTransport.cpp | 29 - Cantera/src/transport/LiquidTransport.h | 22 +- Cantera/src/transport/LiquidTransportData.h | 80 +++ Cantera/src/transport/LiquidTransportParams.h | 3 + Cantera/src/transport/Makefile.in | 6 +- Cantera/src/transport/MixTransport.h | 0 Cantera/src/transport/MultiTransport.h | 0 Cantera/src/transport/SimpleTransport.cpp | 569 +++++++++++++++ Cantera/src/transport/SimpleTransport.h | 666 ++++++++++++++++++ Cantera/src/transport/TransportBase.h | 39 +- Cantera/src/transport/TransportFactory.cpp | 280 ++++++-- Cantera/src/transport/TransportFactory.h | 44 +- 12 files changed, 1597 insertions(+), 141 deletions(-) create mode 100644 Cantera/src/transport/LiquidTransportData.h mode change 100755 => 100644 Cantera/src/transport/MixTransport.h mode change 100755 => 100644 Cantera/src/transport/MultiTransport.h create mode 100644 Cantera/src/transport/SimpleTransport.cpp create mode 100644 Cantera/src/transport/SimpleTransport.h mode change 100755 => 100644 Cantera/src/transport/TransportBase.h mode change 100755 => 100644 Cantera/src/transport/TransportFactory.cpp diff --git a/Cantera/src/transport/LiquidTransport.cpp b/Cantera/src/transport/LiquidTransport.cpp index 8175c32eb..08b566dea 100644 --- a/Cantera/src/transport/LiquidTransport.cpp +++ b/Cantera/src/transport/LiquidTransport.cpp @@ -443,35 +443,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]. * diff --git a/Cantera/src/transport/LiquidTransport.h b/Cantera/src/transport/LiquidTransport.h index 852641192..da8a1751b 100644 --- a/Cantera/src/transport/LiquidTransport.h +++ b/Cantera/src/transport/LiquidTransport.h @@ -326,19 +326,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). @@ -719,14 +706,15 @@ namespace Cantera { int m_nDim; private: - - /** - * Throw an exception if this method is invoked. + + //! 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..f5b1bd248 --- /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 $ + * + * Copyright 2001 California Institute of Technology + * + */ + +#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 08142b745..9b0d98c94 100644 --- a/Cantera/src/transport/LiquidTransportParams.h +++ b/Cantera/src/transport/LiquidTransportParams.h @@ -6,6 +6,7 @@ #include "ct_defs.h" #include "TransportBase.h" #include "TransportParams.h" +#include "LiquidTransportData.h" #include "xml.h" #include "XML_Writer.h" @@ -81,6 +82,8 @@ namespace Cantera { vector_fp B_k_cond; + std::vector LTData; + }; } 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.h b/Cantera/src/transport/MixTransport.h old mode 100755 new mode 100644 diff --git a/Cantera/src/transport/MultiTransport.h b/Cantera/src/transport/MultiTransport.h old mode 100755 new mode 100644 diff --git a/Cantera/src/transport/SimpleTransport.cpp b/Cantera/src/transport/SimpleTransport.cpp new file mode 100644 index 000000000..ed119bfe5 --- /dev/null +++ b/Cantera/src/transport/SimpleTransport.cpp @@ -0,0 +1,569 @@ +/** + * @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), + 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) + { + } + //================================================================================================ + SimpleTransport::SimpleTransport(const SimpleTransport &right) : + Transport(), + m_nsp(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) + { + /* + * 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; + m_tmin = right.m_tmin; + m_tmax = right.m_tmax; + m_mw = right.m_mw; + + m_Grad_X = right.m_Grad_X; + m_Grad_T = right.m_Grad_T; + m_Grad_V = right.m_Grad_V; + + m_viscSpecies = right.m_viscSpecies; + m_condSpecies = right.m_condSpecies; + m_iStateMF = -1; + m_molefracs = right.m_molefracs; + m_concentrations = right.m_concentrations; + 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) { + + // 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()); + + //save logarithm of pre-exponential for easier computation + + //m_diffcoeffs = tr.diffcoeffs; + + + m_viscSpecies.resize(m_nsp); + m_condSpecies.resize(m_nsp); + + + m_molefracs.resize(m_nsp); + 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_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; + } + } + } + //================================================================================================ + void SimpleTransport::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); + } + } + //================================================================================================ + void SimpleTransport::set_Grad_V(const doublereal* const grad_V) { + for (int a = 0; a < m_nDim; a++) { + m_Grad_V[a] = grad_V[a]; + } + } + //================================================================================================ + 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 = - \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 + */ + 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 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*ldf + k] = -rhon * mw[k] * m_spwork[k] * m_Grad_X[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 (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_viscSpecies[k] = coeff[0] * pow(m_temp,coeff[1]) * exp(-coeff[2]/m_temp); + } + } + 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..45f828564 --- /dev/null +++ b/Cantera/src/transport/SimpleTransport.h @@ -0,0 +1,666 @@ +/** + * + * @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 double Coeff_T_ [4]; + + + //! 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); + + + 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 Mobilities + /*! + * @param mobil + */ + virtual void getMobilities(doublereal* const mobil); + + //! 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_; + + //! 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. + vector m_coeffVisc_Ns; + + //! Pure species thermal conductivities in Arrhenius temperature-dependent form. + /*! + * + */ + vector m_coeffLambda_Ns; + + + //! Pure species viscosities in Arrhenius temperature-dependent form. + vector m_coeffDiff_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_; + + + + doublereal meanMolecularWeight_; + 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/TransportBase.h b/Cantera/src/transport/TransportBase.h old mode 100755 new mode 100644 index dd0d7aae3..02b3d2c90 --- a/Cantera/src/transport/TransportBase.h +++ b/Cantera/src/transport/TransportBase.h @@ -45,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; @@ -159,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 @@ -169,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"); } @@ -273,9 +284,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, @@ -298,9 +309,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"); } @@ -348,20 +363,20 @@ 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("init"); return false; } + { err("initGas"); return false; } /** * Called by TransportFactory to set parameters. */ virtual bool initLiquid( LiquidTransportParams& tr ) - { err("init"); return false; } + { 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 2beffa8e9..bb3dfbcd9 --- 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" @@ -96,6 +97,21 @@ namespace Cantera { }; + /** + * 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 ////////////// @@ -239,6 +255,7 @@ namespace Cantera { 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; @@ -319,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; @@ -342,7 +364,7 @@ namespace Cantera { -/** + /** * Prepare to build a new kinetic-theory-based transport manager * for low-density gases. Uses polynomial fits to Monchick & Mason * collision integrals. @@ -473,13 +495,13 @@ 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) { + const std::vector &transport_database, + thermo_t* thermo, int log_level, LiquidTransportParams& trParam) { // constant mixture attributes trParam.thermo = thermo; @@ -506,7 +528,7 @@ namespace Cantera { XML_Node root, log; getLiquidTransportData(transport_database, log, - trParam.thermo->speciesNames(), trParam); + trParam.thermo->speciesNames(), trParam); //int i, j; //for (i = 0; i < nsp; i++) trParam.poly[i].resize(nsp); @@ -617,7 +639,7 @@ namespace Cantera { // set up Monchick and Mason collision integrals setupMM(flog, transport_database, thermo, mode, log_level, trParam); // do model-specific initialization - tran->init(trParam); + tran->initGas(trParam); #ifdef DEBUG_MODE if (m_verbose) { trParam.xml->XML_close(flog, "transport"); @@ -633,8 +655,8 @@ namespace Cantera { * class and calls setupLiquidTransport(). */ void TransportFactory::initLiquidTransport(Transport* tran, - thermo_t* thermo, - int log_level) { + thermo_t* thermo, + int log_level) { const std::vector & transport_database = thermo->speciesData(); @@ -651,7 +673,7 @@ namespace Cantera { #endif setupLiquidTransport(flog, transport_database, thermo, log_level, trParam); // do model-specific initialization - tran->init(trParam); + tran->initLiquid(trParam); #ifdef DEBUG_MODE if (m_verbose) { trParam.xml->XML_close(flog, "transport"); @@ -780,7 +802,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; + // 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. @@ -883,13 +905,17 @@ namespace Cantera { const std::vector &names, LiquidTransportParams& trParam) { - string name; + 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; + 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 @@ -900,53 +926,189 @@ 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; + 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 { - XML_Node& trNode = sp.child("transport"); + if (sp.hasChild("transport")) { + XML_Node& trNode = sp.child("transport"); - hydrodynamic_radius = getFloat(trNode, "hydrodynamic_radius"); + // Fill datatable with LiquidTransportData objects for error checking + // and then insertion into LiquidTransportData objects below. + LiquidTransportData data; + data.speciesName = name; - XML_Node& visc = trNode.child("viscosity"); - getArrhenius(visc, A_visc, n_visc, Tact_visc ); + /* + * 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); + } + } - XML_Node& thermCond = trNode.child("thermal_conductivity"); - getArrhenius(thermCond, A_thcond, n_thcond, Tact_thcond ); + /* + * 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"]); + } + } - // Fill datatable with LiquidTransportData objects for error checking - // and then insertion into LiquidTransportData objects below. - LiquidTransportData data; - data.speciesName = name; + /* + * 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"]); + } + } - if ( hydrodynamic_radius > 0.0) data.hydroradius = hydrodynamic_radius; - else throw TransportDBError(linenum, - "negative or zero hydrodynamic radius"); + + /* + * 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"]); + } + } - if (A_visc >= 0.0) { - data.viscCoeffs[0] = A_visc; - data.viscCoeffs[1] = n_visc; - data.viscCoeffs[2] = Tact_visc; + datatable[name] = data; } - else throw TransportDBError(linenum, - "negative pre-exponential for viscosity"); - - if (A_thcond >= 0.0) { - data.thermalCondCoeffs[0] = A_thcond; - data.thermalCondCoeffs[1] = n_thcond; - data.thermalCondCoeffs[2] = Tact_thcond; - } - else throw TransportDBError(linenum, - "negative pre-exponential for thermalCondoctivity"); - - datatable[name] = data; } catch(CanteraError) { ; } } + trParam.LTData.clear(); for (i = 0; i < trParam.nsp_; i++) { LiquidTransportData& trdat = datatable[names[i]]; @@ -960,24 +1122,36 @@ namespace Cantera { } // parameters should be converted to SI units before storing - - trParam.visc_A[i] = trdat.viscCoeffs[0] ; - trParam.visc_n[i] = trdat.viscCoeffs[1] ; - trParam.visc_Tact[i] = trdat.viscCoeffs[2] ; + 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] ; + } - trParam.thermCond_A[i] = trdat.thermalCondCoeffs[0] ; - trParam.thermCond_n[i] = trdat.thermalCondCoeffs[1] ; - trParam.thermCond_Tact[i] = trdat.thermalCondCoeffs[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; + // 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; } diff --git a/Cantera/src/transport/TransportFactory.h b/Cantera/src/transport/TransportFactory.h index 489ef7156..c89f1b7d1 100755 --- 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,20 +62,12 @@ namespace Cantera { doublereal rotRelaxNumber; }; - struct LiquidTransportData { - LiquidTransportData() : speciesName("-"), - hydroradius(-1) {} - std::string speciesName; - doublereal hydroradius; - vector_fp viscCoeffs; - vector_fp thermalCondCoeffs; - }; // forward references class MMCollisionInt; class GasTransportParams; class LiquidTransportParams; - class XML_Node; + class XML_Node; //! The purpose of TransportFactory is to create new instances of @@ -173,9 +166,19 @@ namespace Cantera { XML_Node& log, const std::vector& names, 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); + XML_Node& log, const std::vector& names, + LiquidTransportParams& tr); /** Generate polynomial fits to viscosity, conductivity, and * binary diffusion coefficients */ @@ -193,8 +196,8 @@ namespace Cantera { void setupLiquidTransport(std::ostream &flog, const std::vector &transport_database, - thermo_t* thermo, int log_level, - LiquidTransportParams& tr); + thermo_t* thermo, int log_level, + LiquidTransportParams& tr); /// Second-order correction to the binary diffusion coefficients @@ -208,21 +211,6 @@ namespace Cantera { const GasTransportParams& tr, doublereal& f_eps, doublereal& f_sigma); - /** - * 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; - } - //! Boolean indicating whether to turn on verbose printing bool m_verbose; From 787fd3778e142a444e0f4062718a1d7f9bbd85ad Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Sat, 10 Oct 2009 21:54:55 +0000 Subject: [PATCH 19/24] More work on the SimpleTransport object. Nominally, the object is finished. But, it needs testing. --- Cantera/src/transport/AqueousTransport.cpp | 20 +- Cantera/src/transport/AqueousTransport.h | 39 ++- Cantera/src/transport/LiquidTransport.cpp | 77 ++++-- Cantera/src/transport/LiquidTransport.h | 37 ++- Cantera/src/transport/LiquidTransportData.h | 4 +- Cantera/src/transport/SimpleTransport.cpp | 256 +++++++++++++++++--- Cantera/src/transport/SimpleTransport.h | 72 +++++- Cantera/src/transport/TransportBase.h | 56 ++++- Cantera/src/transport/TransportFactory.h | 0 9 files changed, 466 insertions(+), 95 deletions(-) mode change 100755 => 100644 Cantera/src/transport/TransportFactory.h diff --git a/Cantera/src/transport/AqueousTransport.cpp b/Cantera/src/transport/AqueousTransport.cpp index 171b256a6..9fda00393 100644 --- a/Cantera/src/transport/AqueousTransport.cpp +++ b/Cantera/src/transport/AqueousTransport.cpp @@ -181,6 +181,7 @@ namespace Cantera { /******************* binary diffusion coefficients **************/ + //================================================================================================ void AqueousTransport::getBinaryDiffCoeffs(const int ld, doublereal* const d) { int i,j; @@ -197,32 +198,27 @@ namespace Cantera { d[ld*j + i] = rp * m_bdiff(i,j); } } - - + //================================================================================================ 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); + 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++) { diff --git a/Cantera/src/transport/AqueousTransport.h b/Cantera/src/transport/AqueousTransport.h index 07e3e3dd9..1d593bd16 100644 --- a/Cantera/src/transport/AqueousTransport.h +++ b/Cantera/src/transport/AqueousTransport.h @@ -194,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 /*! diff --git a/Cantera/src/transport/LiquidTransport.cpp b/Cantera/src/transport/LiquidTransport.cpp index 08b566dea..903b1298f 100644 --- a/Cantera/src/transport/LiquidTransport.cpp +++ b/Cantera/src/transport/LiquidTransport.cpp @@ -305,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 **********************/ /* @@ -652,7 +690,7 @@ namespace Cantera { */ void LiquidTransport::updateCond_temp() { - int k; + /* if (m_mode == CK_Mode) { for (k = 0; k < m_nsp; k++) { @@ -676,8 +714,7 @@ 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++) { diff --git a/Cantera/src/transport/LiquidTransport.h b/Cantera/src/transport/LiquidTransport.h index da8a1751b..c3df40c3f 100644 --- a/Cantera/src/transport/LiquidTransport.h +++ b/Cantera/src/transport/LiquidTransport.h @@ -251,11 +251,42 @@ namespace Cantera { */ virtual doublereal thermalConductivity(); - //! Get the Mobilities + //! Get the Electrical mobilities (m^2/V/s). /*! - * @param mobil + * 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); + 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 /*! diff --git a/Cantera/src/transport/LiquidTransportData.h b/Cantera/src/transport/LiquidTransportData.h index f5b1bd248..a664bedc1 100644 --- a/Cantera/src/transport/LiquidTransportData.h +++ b/Cantera/src/transport/LiquidTransportData.h @@ -8,7 +8,7 @@ * $Date: 2008/12/24 18:19:01 $ * $Revision: 1.14 $ * - * Copyright 2001 California Institute of Technology + * * */ @@ -63,7 +63,7 @@ namespace Cantera { //! Model type for the hydroradius LiquidTR_Model model_viscosity; - vector_fp viscCoeffs; + vector_fp viscCoeffs; //! Model type for the hydroradius LiquidTR_Model model_thermalCond; diff --git a/Cantera/src/transport/SimpleTransport.cpp b/Cantera/src/transport/SimpleTransport.cpp index ed119bfe5..618f6fc3a 100644 --- a/Cantera/src/transport/SimpleTransport.cpp +++ b/Cantera/src/transport/SimpleTransport.cpp @@ -31,9 +31,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), @@ -43,13 +48,18 @@ namespace Cantera { m_diff_mix_ok(false), m_diff_temp_ok(false), m_cond_temp_ok(false), - m_cond_mix_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), @@ -62,7 +72,8 @@ namespace Cantera { m_diff_mix_ok(false), m_diff_temp_ok(false), m_cond_temp_ok(false), - m_cond_mix_ok(false) + m_cond_mix_ok(false), + m_nDim(1) { /* * Use the assignment operator to do the brunt @@ -76,23 +87,36 @@ namespace Cantera { 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; @@ -120,7 +144,7 @@ namespace Cantera { * 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(); @@ -132,21 +156,122 @@ namespace Cantera { copy(m_thermo->molecularWeights().begin(), m_thermo->molecularWeights().end(), m_mw.begin()); - //save logarithm of pre-exponential for easier computation - - //m_diffcoeffs = tr.diffcoeffs; - - + /* + * Get the input Viscosities + */ m_viscSpecies.resize(m_nsp); - m_condSpecies.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); @@ -232,19 +357,64 @@ namespace Cantera { } } //================================================================================================ + // 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) { - // this needs to be checked out. int k; getMixDiffCoeffs(DATA_PTR(m_spwork)); - doublereal c1 = ElectronCharge / (Boltzmann * m_temp); + doublereal t = m_thermo->temperature(); + doublereal c1 = ElectronCharge / (Boltzmann * t); 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 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; } } //================================================================================================ @@ -260,7 +430,6 @@ namespace Cantera { m_Grad_X[i] = grad_X[i]; } } - //================================================================================================ // Returns the mixture thermal conductivity of the solution /* @@ -333,7 +502,7 @@ namespace Cantera { } //================================================================================================ // Return the species diffusive mass fluxes wrt to - // the mass averaged velocity, + // the mass averaged velocity. /* * * units = kg/m2/s @@ -346,13 +515,13 @@ namespace Cantera { * formula * * \f[ - * j_k = - \rho M_k D_k \nabla X_k - Y_k V_c + * 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 {\rho M_j D_j \nabla X_j} + * 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 @@ -367,16 +536,28 @@ namespace Cantera { getMixDiffCoeffs(DATA_PTR(m_spwork)); - const array_fp& mw = m_thermo->molecularWeights(); const doublereal* y = m_thermo->massFractions(); - doublereal rhon = m_thermo->molarDensity(); + doublereal conc = 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*ldf + k] = -rhon * mw[k] * m_spwork[k] * m_Grad_X[n*m_nsp + k]; - sum[n] += fluxes[n*ldf + k]; + + 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 @@ -478,22 +659,31 @@ namespace Cantera { */ void SimpleTransport::updateDiff_T() { int k; - if (tempDepType_ == 0) { - for (k = 0; k < m_nsp; k++) { - Coeff_T_ &coeff = m_coeffDiff_Ns[k]; - m_diffSpecies[k] = coeff[0]; + 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 if (tempDepType_ == 1) { + } else { + double visc = viscosity(); + double RT = GasConstant * m_temp; for (k = 0; k < m_nsp; k++) { - Coeff_T_ &coeff = m_coeffDiff_Ns[k]; - m_viscSpecies[k] = coeff[0] * pow(m_temp,coeff[1]) * exp(-coeff[2]/m_temp); + 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. */ @@ -550,7 +740,7 @@ namespace Cantera { return true; } - + //================================================================================================ /** * Throw an exception if this method is invoked. * This probably indicates something is not yet implemented. diff --git a/Cantera/src/transport/SimpleTransport.h b/Cantera/src/transport/SimpleTransport.h index 45f828564..35eb9560f 100644 --- a/Cantera/src/transport/SimpleTransport.h +++ b/Cantera/src/transport/SimpleTransport.h @@ -129,7 +129,7 @@ namespace Cantera { class SimpleTransport : public Transport { public: - typedef double Coeff_T_ [4]; + typedef vector_fp Coeff_T_; //! Default constructor. @@ -237,6 +237,12 @@ namespace Cantera { 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); @@ -252,6 +258,7 @@ namespace Cantera { * Solvent-only: * \f[ * \lambda = \lambda_0 + * \f] * Mixture-average: * \f[ @@ -262,13 +269,46 @@ namespace Cantera { * * @see updateCond_T(); */ + virtual doublereal thermalConductivity(); - //! Get the Mobilities + //! Get the electrical Mobilities (m^2/V/s). /*! - * @param mobil + * 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); + + //! 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 /*! @@ -428,6 +468,15 @@ namespace Cantera { */ 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; @@ -441,17 +490,20 @@ namespace Cantera { vector_fp m_mw; //! Pure species viscosities in Arrhenius temperature-dependent form. - vector m_coeffVisc_Ns; + std::vector m_coeffVisc_Ns; //! Pure species thermal conductivities in Arrhenius temperature-dependent form. /*! * */ - vector m_coeffLambda_Ns; + std::vector m_coeffLambda_Ns; //! Pure species viscosities in Arrhenius temperature-dependent form. - vector m_coeffDiff_Ns; + std::vector m_coeffDiff_Ns; + + + std::vector m_coeffHydroRadius_Ns; //! Internal value of the gradient of the mole fraction vector @@ -577,9 +629,10 @@ namespace Cantera { */ doublereal concTot_; - - + //! Mean molecular weight doublereal meanMolecularWeight_; + + //! Density doublereal dens_; //! Local copy of the charge of each species @@ -587,7 +640,6 @@ namespace Cantera { * Contains the charge of each species (length m_nsp) */ vector_fp m_chargeSpecies; - //! Current Temperature -> locally storred /*! diff --git a/Cantera/src/transport/TransportBase.h b/Cantera/src/transport/TransportBase.h index 02b3d2c90..864469f8c 100644 --- a/Cantera/src/transport/TransportBase.h +++ b/Cantera/src/transport/TransportBase.h @@ -193,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"); } + //@} @@ -258,12 +292,12 @@ namespace Cantera { * 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) { + 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 ); } diff --git a/Cantera/src/transport/TransportFactory.h b/Cantera/src/transport/TransportFactory.h old mode 100755 new mode 100644 From 4ee4b4891fb9a9bb6e131b902f8958eac21a8c92 Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Sun, 11 Oct 2009 19:20:11 +0000 Subject: [PATCH 20/24] Added a missing function to AqueousTransport --- Cantera/src/transport/AqueousTransport.cpp | 25 ++++++++++++++++++++++ Cantera/src/transport/SimpleTransport.cpp | 3 +-- 2 files changed, 26 insertions(+), 2 deletions(-) diff --git a/Cantera/src/transport/AqueousTransport.cpp b/Cantera/src/transport/AqueousTransport.cpp index 9fda00393..55e820462 100644 --- a/Cantera/src/transport/AqueousTransport.cpp +++ b/Cantera/src/transport/AqueousTransport.cpp @@ -199,7 +199,32 @@ namespace Cantera { } } //================================================================================================ + // 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) { + 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]; + } + } + //================================================================================================ + 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++) { diff --git a/Cantera/src/transport/SimpleTransport.cpp b/Cantera/src/transport/SimpleTransport.cpp index 618f6fc3a..eafaa02fe 100644 --- a/Cantera/src/transport/SimpleTransport.cpp +++ b/Cantera/src/transport/SimpleTransport.cpp @@ -376,8 +376,7 @@ namespace Cantera { void SimpleTransport::getMobilities(doublereal* const mobil) { int k; getMixDiffCoeffs(DATA_PTR(m_spwork)); - doublereal t = m_thermo->temperature(); - doublereal c1 = ElectronCharge / (Boltzmann * t); + doublereal c1 = ElectronCharge / (Boltzmann * m_temp); for (k = 0; k < m_nsp; k++) { mobil[k] = c1 * m_spwork[k]; } From 5fd13d3591ca984dc7a1093f49522d8a5a33099a Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Sun, 11 Oct 2009 20:46:57 +0000 Subject: [PATCH 21/24] Took out default parameters for some of the transport initialization routines. It makes no sense to call these things if you don't have a ThermoPhase already created. --- Cantera/src/transport/TransportFactory.cpp | 2 +- Cantera/src/transport/TransportFactory.h | 6 +++--- 2 files changed, 4 insertions(+), 4 deletions(-) diff --git a/Cantera/src/transport/TransportFactory.cpp b/Cantera/src/transport/TransportFactory.cpp index bb3dfbcd9..50c324461 100644 --- a/Cantera/src/transport/TransportFactory.cpp +++ b/Cantera/src/transport/TransportFactory.cpp @@ -356,7 +356,7 @@ namespace Cantera { break; #endif default: - throw CanteraError("newTransport","unknown transport model"); + throw CanteraError("newTransport","unknown transport model: " + transportModel); } phase->restoreState(state); return tr; diff --git a/Cantera/src/transport/TransportFactory.h b/Cantera/src/transport/TransportFactory.h index c89f1b7d1..299e4397f 100644 --- a/Cantera/src/transport/TransportFactory.h +++ b/Cantera/src/transport/TransportFactory.h @@ -128,15 +128,15 @@ namespace Cantera { /// Build a new transport manager 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); /// 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); From d503c55864d3233eec7206215c9a93e84457d946 Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Sun, 11 Oct 2009 21:34:06 +0000 Subject: [PATCH 22/24] Added a clear() function to the xml object. This will wipe out the current level and everything below it. --- Cantera/src/base/xml.cpp | 21 +++++++++++++++++++++ Cantera/src/base/xml.h | 9 ++++++++- 2 files changed, 29 insertions(+), 1 deletion(-) 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 From a258cdf0839c598726c81d9f0aeb484d282788c5 Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Sun, 11 Oct 2009 21:38:42 +0000 Subject: [PATCH 23/24] Added an xml storage section to importPhase(). Note this represents a long standing bug in some ThermoPhase implementations, where the xml tree was not being storred for the phase node in some cases. --- Cantera/src/thermo/ThermoFactory.cpp | 11 +++++++++++ 1 file changed, 11 insertions(+) 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()); From 3b4401a77d73d4f3ce760bc08d49aeac387ae2ff Mon Sep 17 00:00:00 2001 From: Harry Moffat Date: Mon, 12 Oct 2009 14:13:10 +0000 Subject: [PATCH 24/24] Added into Transport Factory an option to read the phase XML file and implement the transport manager listed there. --- Cantera/src/numerics/Integrator.h | 13 +++++- Cantera/src/transport/TransportFactory.cpp | 23 ++++++++++- Cantera/src/transport/TransportFactory.h | 46 +++++++++++++++++----- 3 files changed, 70 insertions(+), 12 deletions(-) diff --git a/Cantera/src/numerics/Integrator.h b/Cantera/src/numerics/Integrator.h index cd442db88..2006f62e8 100644 --- 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,7 +83,9 @@ 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. /*! diff --git a/Cantera/src/transport/TransportFactory.cpp b/Cantera/src/transport/TransportFactory.cpp index 50c324461..3f409b8db 100644 --- a/Cantera/src/transport/TransportFactory.cpp +++ b/Cantera/src/transport/TransportFactory.cpp @@ -295,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; @@ -362,6 +362,27 @@ namespace Cantera { 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); + } /** diff --git a/Cantera/src/transport/TransportFactory.h b/Cantera/src/transport/TransportFactory.h index 299e4397f..d38603557 100644 --- a/Cantera/src/transport/TransportFactory.h +++ b/Cantera/src/transport/TransportFactory.h @@ -124,12 +124,26 @@ 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, 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, int mode=0, int log_level=0); @@ -139,8 +153,6 @@ namespace Cantera { thermo_t* thermo, int log_level=0); - - private: //! Static instance of the factor -> This is the only instance of this @@ -150,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 /*! @@ -228,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(); @@ -241,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; }