From f312f0d1d0cf52f9d2c9827a9f0a9c4b3b0d1abb Mon Sep 17 00:00:00 2001 From: ignis Date: Fri, 6 Mar 2020 11:12:43 +0900 Subject: [PATCH] WmeanTransport: Transport properties are calculated as a weighted mean of properties of selected species --- include/cantera/transport/WmeanTransport.h | 106 ++++++++++++ src/transport/TransportFactory.cpp | 3 + src/transport/WmeanTransport.cpp | 177 +++++++++++++++++++++ 3 files changed, 286 insertions(+) create mode 100644 include/cantera/transport/WmeanTransport.h create mode 100644 src/transport/WmeanTransport.cpp diff --git a/include/cantera/transport/WmeanTransport.h b/include/cantera/transport/WmeanTransport.h new file mode 100644 index 000000000..2a4b80324 --- /dev/null +++ b/include/cantera/transport/WmeanTransport.h @@ -0,0 +1,106 @@ +/** + * @file WmeanTransport.h + * Headers for the WmeanTransport object, which models transport + * properties in ideal gas solutions using the unity Lewis number + * approximation + * (see \ref tranprops and \link Cantera::WmeanTransport WmeanTransport \endlink) . + */ + +// This file is part of Cantera. See License.txt in the top-level directory or +// at https://cantera.org/license.txt for license and copyright information. + +#ifndef CT_WMEANTRAN_H +#define CT_WMEANTRAN_H + +#include "MixTransport.h" + +namespace Cantera +{ +//! Class WmeanTransport implements the unity Lewis number approximation +//! for the mixture-averaged species diffusion coefficients. Mixture-averaged +//! transport properties for viscosity and thermal conductivity are inherited +//! from the MixTransport class. +//! @ingroup tranprops +class WmeanTransport : public MixTransport +{ +public: + WmeanTransport(); + + virtual std::string transportType() const { + return "Wmean"; + } + + virtual void init(thermo_t* thermo, int mode=0, int log_level=0); + + //! Returns the unity Lewis number approximation based diffusion + //! coefficients [m^2/s]. + /*! + * Returns the unity Lewis number approximation based diffusion coefficients + * for a gas, appropriate for calculating the mass averaged diffusive flux + * with respect to the mass averaged velocity using gradients of the mole + * fraction. + * + * \f[ + * D^\prime_{km} = \frac{\lambda}{\rho c_p} + * \f] + * + * In order to obtain the expected behavior from a unity Lewis number model, + * this formulation requires that the correction velocity be computed as + * + * \f[ + * V_c = \sum \frac{W_k}{\overline{W}} D^\prime_{km} \nabla X_k + * \f] + * + * @param[out] d Vector of diffusion coefficients for each species (m^2/s). + * length m_nsp. + */ + virtual void getMixDiffCoeffs(double* const d) { + GasTransport::getMixDiffCoeffs(d); + + double Dm = 0.0; + for (size_t k = 0; k < m_nsp; k++) { + Dm += d[k] * m_weight[k]; + } + for (size_t k = 0; k < m_nsp; k++) { + d[k] = Dm; + } + } + + //! Not implemented for unity Lewis number approximation + virtual void getMixDiffCoeffsMole(double* const d){ + throw NotImplementedError("WmeanTransport::getMixDiffCoeffsMole"); + } + + //! Returns the unity Lewis number approximation based diffusion + //! coefficients [m^2/s]. + /*! + * These are the coefficients for calculating the diffusive mass fluxes + * from the species mass fraction gradients, computed as + * + * \f[ + * D_{km} = \frac{\lambda}{\rho c_p} + * \f] + * + * @param[out] d Vector of diffusion coefficients for each species (m^2/s). + * length m_nsp. + */ + virtual void getMixDiffCoeffsMass(double* const d){ + GasTransport::getMixDiffCoeffsMass(d); + + double Dm = 0.0; + for (size_t k = 0; k < m_nsp; k++) { + Dm += d[k] * m_weight[k]; + } + for (size_t k = 0; k < m_nsp; k++) { + d[k] = Dm; + } + } + +protected: + + //! Internal storage for weights to calculate mean diffusivity + vector_fp m_weight; + +}; +} +#endif diff --git a/src/transport/TransportFactory.cpp b/src/transport/TransportFactory.cpp index 146f69f86..58d071f6d 100644 --- a/src/transport/TransportFactory.cpp +++ b/src/transport/TransportFactory.cpp @@ -8,6 +8,7 @@ #include "cantera/transport/MixTransport.h" #include "cantera/transport/UnityLewisTransport.h" #include "cantera/transport/ConstantTransport.h" +#include "cantera/transport/WmeanTransport.h" #include "cantera/transport/IonGasTransport.h" #include "cantera/transport/WaterTransport.h" #include "cantera/transport/DustyGasTransport.h" @@ -48,6 +49,8 @@ TransportFactory::TransportFactory() m_synonyms["None"] = ""; reg("Constant", []() { return new ConstantTransport(); }); m_synonyms["constant"] = "Constant"; + reg("Wmean", []() { return new WmeanTransport(); }); + m_synonyms["wmean"] = "Wmean"; reg("UnityLewis", []() { return new UnityLewisTransport(); }); m_synonyms["unity-Lewis-number"] = "UnityLewis"; reg("Mix", []() { return new MixTransport(); }); diff --git a/src/transport/WmeanTransport.cpp b/src/transport/WmeanTransport.cpp new file mode 100644 index 000000000..d128e2ff2 --- /dev/null +++ b/src/transport/WmeanTransport.cpp @@ -0,0 +1,177 @@ +/** + * @file WmeanTransport.cpp + * Mixture-averaged transport properties for ideal gas mixtures. + */ + +// This file is part of Cantera. See License.txt in the top-level directory or +// at https://cantera.org/license.txt for license and copyright information. + +#include "cantera/transport/WmeanTransport.h" +#include "cantera/base/stringUtils.h" + +using namespace std; + +namespace Cantera +{ +WmeanTransport::WmeanTransport() : + MixTransport(), + m_weight(m_nsp) +{ +} + +void WmeanTransport::init(ThermoPhase* thermo, int mode, int log_level) +{ + MixTransport::init(thermo, mode, log_level); + + m_weight.resize(m_nsp); + + m_weight.assign(m_nsp, 0.0); + m_weight[0] = 1.0; + + ifstream myfile; + myfile.open ("wmean-weights.txt"); + for (size_t k = 0; k < m_nsp; k++) { + myfile >> m_weight[k]; + + cout << m_thermo->speciesName(k) << m_weight[k] << endl; + } + myfile.close(); + + // m_reftemp = thermo->temperature(); + // cout << "check reference temperature" << m_reftemp << endl; +} + + /* +void WmeanTransport::init(ThermoPhase* thermo, int mode, int log_level) +{ + MixTransport::init(thermo, mode, log_level); +} + +void WmeanTransport::getMobilities(doublereal* const mobil) +{ + getMixDiffCoeffs(m_spwork.data()); + doublereal c1 = ElectronCharge / (Boltzmann * m_temp); + for (size_t k = 0; k < m_nsp; k++) { + mobil[k] = c1 * m_spwork[k]; + } +} + +doublereal WmeanTransport::thermalConductivity() +{ + update_T(); + update_C(); + if (!m_spcond_ok) { + updateCond_T(); + } + if (!m_condmix_ok) { + doublereal sum1 = 0.0, sum2 = 0.0; + for (size_t k = 0; k < m_nsp; k++) { + sum1 += m_molefracs[k] * m_cond[k]; + sum2 += m_molefracs[k] / m_cond[k]; + } + m_lambda = 0.5*(sum1 + 1.0/sum2); + m_condmix_ok = true; + } + return m_lambda; +} + +void WmeanTransport::getThermalDiffCoeffs(doublereal* const dt) +{ + for (size_t k = 0; k < m_nsp; k++) { + dt[k] = 0.0; + } +} + +void WmeanTransport::getSpeciesFluxes(size_t ndim, const doublereal* const grad_T, + size_t ldx, const doublereal* const grad_X, + size_t ldf, doublereal* const fluxes) +{ + update_T(); + update_C(); + getMixDiffCoeffs(m_spwork.data()); + const vector_fp& mw = m_thermo->molecularWeights(); + const doublereal* y = m_thermo->massFractions(); + doublereal rhon = m_thermo->molarDensity(); + vector_fp sum(ndim,0.0); + for (size_t n = 0; n < ndim; n++) { + for (size_t k = 0; k < m_nsp; k++) { + fluxes[n*ldf + k] = -rhon * mw[k] * m_spwork[k] * grad_X[n*ldx + k]; + sum[n] += fluxes[n*ldf + k]; + } + } + // add correction flux to enforce sum to zero + for (size_t n = 0; n < ndim; n++) { + for (size_t k = 0; k < m_nsp; k++) { + fluxes[n*ldf + k] -= y[k]*sum[n]; + } + } +} + +void WmeanTransport::update_T() +{ + doublereal t = m_thermo->temperature(); + if (t == m_temp && m_nsp == m_thermo->nSpecies()) { + return; + } + if (t < 0.0) { + throw CanteraError("WmeanTransport::update_T", + "negative temperature {}", t); + } + GasTransport::update_T(); + // temperature has changed, so polynomial fits will need to be redone. + m_spcond_ok = false; + m_bindiff_ok = false; + m_condmix_ok = false; +} + +void WmeanTransport::update_C() +{ + // signal that concentration-dependent quantities will need to be recomputed + // before use, and update the local mole fractions. + m_visc_ok = false; + m_condmix_ok = false; + m_thermo->getMoleFractions(m_molefracs.data()); + + // add an offset to avoid a pure species condition + for (size_t k = 0; k < m_nsp; k++) { + m_molefracs[k] = std::max(Tiny, m_molefracs[k]); + } +} + +void WmeanTransport::updateCond_T() +{ + if (m_mode == CK_Mode) { + for (size_t k = 0; k < m_nsp; k++) { + m_cond[k] = exp(dot4(m_polytempvec, m_condcoeffs[k])); + } + } else { + for (size_t k = 0; k < m_nsp; k++) { + m_cond[k] = m_sqrt_t * dot5(m_polytempvec, m_condcoeffs[k]); + } + } + m_spcond_ok = true; + m_condmix_ok = false; +} + +void WmeanTransport::update_T() +{ + doublereal t = m_thermo->temperature(); + if (m_reftemp == m_temp && m_nsp == m_thermo->nSpecies()) { + return; + } + if (t < 0.0) { + throw CanteraError("WmeanTransport::update_T", + "negative temperature {}", t); + } + m_thermo->setTemperature(m_reftemp); + GasTransport::update_T(); + m_thermo->setTemperature(t); + // temperature has changed, so polynomial fits will need to be redone. + m_spcond_ok = false; + m_bindiff_ok = false; + m_condmix_ok = false; +} + + */ + +}