research/cantera
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Got most of the way towards creating a new transport Object SimpleTransport

Browse source 
that obeys the simplest formulas possible.
This commit is contained in:
Harry Moffat 2009-10-10 04:23:04 +00:00
parent 636c63482b
commit c1d59a912b
12 changed files with 1597 additions and 141 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

29
Cantera/src/transport/LiquidTransport.cpp
View file

@ -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].
*

22
Cantera/src/transport/LiquidTransport.h
View file

@ -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

80
Cantera/src/transport/LiquidTransportData.h Normal file
View file

@ -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 <vector>
#include <string>
#include <iostream>
#include <new>
// 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

3
Cantera/src/transport/LiquidTransportParams.h
View file

@ -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<Cantera::LiquidTransportData> LTData;
};
}

6
Cantera/src/transport/Makefile.in
View file

@ -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

0
Cantera/src/transport/MixTransport.h Executable file → Normal file
View file

0
Cantera/src/transport/MultiTransport.h Executable file → Normal file
View file

569
Cantera/src/transport/SimpleTransport.cpp Normal file
View file

@ -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 <iostream>
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<Transport *>(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;
}
//================================================================================================
}
//================================================================================================

666
Cantera/src/transport/SimpleTransport.h Normal file
View file

@ -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 <vector>
#include <string>
#include <map>
#include <numeric>
#include <algorithm>
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.
*
* <H2> Viscosity Calculation </H2>
*
* 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<Coeff_T_> m_coeffVisc_Ns;
//! Pure species thermal conductivities in Arrhenius temperature-dependent form.
/*!
*
*/
vector<Coeff_T_> m_coeffLambda_Ns;
//! Pure species viscosities in Arrhenius temperature-dependent form.
vector<Coeff_T_> 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

39
Cantera/src/transport/TransportBase.h Executable file → Normal file
View file

@ -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; }

280
Cantera/src/transport/TransportFactory.cpp Executable file → Normal file
View file

@ -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<const XML_Node*> &transport_database,
thermo_t* thermo, int log_level, LiquidTransportParams& trParam) {
const std::vector<const XML_Node*> &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<const XML_Node*> & 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<std::string> &names,
LiquidTransportParams& trParam)
{
string name;
std::string name;
/*
* Create a map of species names versus liquid transport data parameters
*/
std::map<std::string, LiquidTransportData> 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<int>(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:
* <hydrodynamic_radius model="Constant"> 3.0 </hydrodynamic_radius>
* <hydrodynamic_radius> 3.0 </hydrodynamic_radius>
*/
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:
* <viscosity model="Constant"> 3.0 </viscosity>
* <viscosity> 3.0 </viscosity>
* <viscosity model="Arrhenius">
* <A units="Pa S"> 1.0 </A>
* <b> 2.0 </b>
* <E units="kcal/gmol"> 3.0 </E>
* </viscosity>
*
* <viscosity model="Coeff">
* <float_array> 0.0. 1.0, 2.0, 3.0, 4.0 </float_array>
* </viscosity>
*
*/
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:
* <thermal_conductivity model="Constant"> 3.0 </thermal_conductivity>
* <thermal_conductivity> 3.0 </thermal_conductivity>
* <thermal_conductivity model="Arrhenius">
* <A units="Pa S"> 1.0 </A>
* <b> 2.0 </b>
* <E units="kcal/gmol"> 3.0 </E>
* </thermal_conductivity>
*
* <thermal_conductivity model="Coeff">
* <float_array> 0.0. 1.0, 2.0, 3.0, 4.0 </float_array>
* </thermal_conductivity>
*
*/
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:
* <speciesDiffusivity model="Constant"> 3.0 </speciesDiffusivity>
* <speciesDiffusivity> 3.0 </speciesDiffusivity>
* <speciesDiffusivity model="Arrhenius">
* <A units="Pa S"> 1.0 </A>
* <b> 2.0 </b>
* <E units="kcal/gmol"> 3.0 </E>
* </speciesDiffusivity>
*
* <speciesDiffusivity model="Coeff">
* <float_array> 0.0. 1.0, 2.0, 3.0, 4.0 </float_array>
* </speciesDiffusivity>
*
*/
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;
}

44
Cantera/src/transport/TransportFactory.h
View file

@ -34,6 +34,7 @@
#include "ct_defs.h"
#include "TransportBase.h"
#include "FactoryBase.h"
#include "LiquidTransportData.h"
#if defined(THREAD_SAFE_CANTERA)
#include <boost/thread/mutex.hpp>
@ -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<std::string>& 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<const XML_Node*> &db,
XML_Node& log, const std::vector<std::string>& names,
LiquidTransportParams& tr);
XML_Node& log, const std::vector<std::string>& 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<const XML_Node*> &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;