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

LiquidTransport.h LiquidTransport.cpp

Browse source 
Changed order of some prototypes in LiquidTransport.h to group public,
protected and private members better.

Renamed methods (consistency with other transport models) like 
virtual void update_temp() 
to
virtual bool update_T()  //note boolean return value

Added method getSpeciesHydrodynamicRadius(doublereal* const radius);
to fill the hydrodynamic radius array.

Added protected member vectors to hold information on
temperature-dependent species transport properties:
(1) model enumeration like  
vector<LiquidTR_Model> m_viscTempDepType_Ns;
(2) coefficients of temperature-dependent properties like
vector<Coeff_T_>  m_coeffVisc_Ns;  
These replace sets of variables like m_visc_A, m_visc_n, m_visc_Tact.
These are filled in LiquidTransport::initLiquid using
LiquidTransportData vector held within LiquidTransportParams

Species thermal conductivity variables now use "lambda" instead of
"cond" for simplicity.

Still need to work on species-species interactions.  Variables like
int m_compositionDepType will need to be replaced by property-specific
(i.e. viscosity, thermal conductivity) variables.  

Major changes to LiquidTransport::initLiquid(LiquidTransportParams&
tr)
- starting to bring in species-species interactions, but not there
yet.  
- complete rewrite of filling temperture-dependent property variables.
- added hydrodynamic radius and species diffusivity property parsing.
- 

LiquidTransport::thermalConductivity() now uses mass-fraction weighted
mixing rule.  Generality will follow.

LiquidTransport::updateViscosity
LiquidTransport::updateCond_T() hava been updated to use new
temperture-dependent property coefficients (m_lambdaTempDepType_Ns and
m_coeffLambda_Ns).

LiquidTransport::updateDiff() now computes Stefan-Maxwell interaction
parameters (D_ij) using Stokes-Einstein reltion.
This commit is contained in:
John Hewson 2009-10-23 18:59:18 +00:00
parent 2492db9750
commit dbef25a481
2 changed files with 556 additions and 214 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

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

@ -36,6 +36,7 @@ namespace Cantera {
m_nsp(0), m_nsp(0),
m_tmin(-1.0), m_tmin(-1.0),
m_tmax(100000.), m_tmax(100000.),
m_compositionDepType(-1),
m_iStateMF(-1), m_iStateMF(-1),
m_temp(-1.0), m_temp(-1.0),
m_logt(0.0), m_logt(0.0),
@ -60,6 +61,7 @@ namespace Cantera {
m_nsp(0), m_nsp(0),
m_tmin(-1.0), m_tmin(-1.0),
m_tmax(100000.), m_tmax(100000.),
m_compositionDepType(-1),
m_iStateMF(-1), m_iStateMF(-1),
m_temp(-1.0), m_temp(-1.0),
m_logt(0.0), m_logt(0.0),
@ -92,17 +94,17 @@ namespace Cantera {
m_tmin = right.m_tmin; m_tmin = right.m_tmin;
m_tmax = right.m_tmax; m_tmax = right.m_tmax;
m_mw = right.m_mw; m_mw = right.m_mw;
m_visc_A = right.m_visc_A; m_viscTempDepType_Ns = right.m_viscTempDepType_Ns;
m_visc_logA = right.m_visc_logA; m_lambdaTempDepType_Ns = right.m_lambdaTempDepType_Ns;
m_visc_n = right.m_visc_n; m_diffTempDepType_Ns = right.m_diffTempDepType_Ns;
m_visc_Tact = right.m_visc_Tact; m_radiusTempDepType_Ns = right.m_radiusTempDepType_Ns;
m_coeffVisc_Ns = right.m_coeffVisc_Ns;
m_coeffLambda_Ns = right.m_coeffLambda_Ns;
m_coeffDiff_Ns = right.m_coeffDiff_Ns;
m_coeffRadius_Ns = right.m_coeffRadius_Ns;
m_visc_Eij = right.m_visc_Eij; m_visc_Eij = right.m_visc_Eij;
m_visc_Sij = right.m_visc_Sij; 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_hydrodynamic_radius = right.m_hydrodynamic_radius;
m_diffcoeffs = right.m_diffcoeffs;
m_Grad_X = right.m_Grad_X; m_Grad_X = right.m_Grad_X;
m_Grad_T = right.m_Grad_T; m_Grad_T = right.m_Grad_T;
m_Grad_V = right.m_Grad_V; m_Grad_V = right.m_Grad_V;
@ -110,7 +112,7 @@ namespace Cantera {
m_bdiff = right.m_bdiff; m_bdiff = right.m_bdiff;
m_viscSpecies = right.m_viscSpecies; m_viscSpecies = right.m_viscSpecies;
m_logViscSpecies = right.m_logViscSpecies; m_logViscSpecies = right.m_logViscSpecies;
m_condSpecies = right.m_condSpecies; m_lambdaSpecies = right.m_lambdaSpecies;
m_iStateMF = -1; m_iStateMF = -1;
m_molefracs = right.m_molefracs; m_molefracs = right.m_molefracs;
m_concentrations = right.m_concentrations; m_concentrations = right.m_concentrations;
@ -134,7 +136,6 @@ namespace Cantera {
m_cond_temp_ok = false; m_cond_temp_ok = false;
m_cond_mix_ok = false; m_cond_mix_ok = false;
m_mode = right.m_mode; m_mode = right.m_mode;
m_diam = right.m_diam;
m_debug = right.m_debug; m_debug = right.m_debug;
m_nDim = right.m_nDim; m_nDim = right.m_nDim;
@ -153,45 +154,214 @@ namespace Cantera {
*/ */
bool LiquidTransport::initLiquid(LiquidTransportParams& tr) { bool LiquidTransport::initLiquid(LiquidTransportParams& tr) {
int k;
// constant substance attributes // constant substance attributes
m_thermo = tr.thermo; m_thermo = tr.thermo;
m_nsp = m_thermo->nSpecies(); m_nsp = m_thermo->nSpecies();
m_tmin = m_thermo->minTemp(); m_tmin = m_thermo->minTemp();
m_tmax = m_thermo->maxTemp(); m_tmax = m_thermo->maxTemp();
/*
* Read the transport block in the phase XML Node
* It's not an error if this block doesn't exist. Just use the defaults
*/
XML_Node &phaseNode = m_thermo->xml();
if (phaseNode.hasChild("transport")) {
XML_Node& transportNode = phaseNode.child("transport");
if ( transportNode.hasChild("viscosity")) {
XML_Node& viscosityNode = transportNode.child("viscosity");
string viscosityModel = viscosityNode.attrib("model");
if (viscosityModel == "") {
throw CanteraError("LiquidTransport::initLiquid",
"transport::visosity XML node doesn't have a model string");
}
}
string transportModel = transportNode.attrib("model");
if (transportModel == "LiquidTransport") {
/*
* <compositionDependence model="Solvent_Only"/>
* or
* <compositionDependence model="Mixture_Averaged"/>
*/
std::string modelName = "";
if (getOptionalModel(transportNode, "compositionDependence",
modelName)) {
modelName = lowercase(modelName);
if (modelName == "solvent_only") {
m_compositionDepType = 0;
} else if (modelName == "mixture_averaged") {
m_compositionDepType = 1;
} else {
throw CanteraError("LiquidTransport::initLiquid", "Unknown compositionDependence Model: " + modelName);
}
}
}
}
// make a local copy of the molecular weights // make a local copy of the molecular weights
m_mw.resize(m_nsp); m_mw.resize(m_nsp);
copy(m_thermo->molecularWeights().begin(), copy(m_thermo->molecularWeights().begin(),
m_thermo->molecularWeights().end(), m_mw.begin()); m_thermo->molecularWeights().end(), m_mw.begin());
// copy parameters into local storage /*
m_visc_A = tr.visc_A ; * Get the input Viscosities
m_visc_n = tr.visc_n ; */
m_visc_Tact = tr.visc_Tact ; m_viscSpecies.resize(m_nsp);
m_coeffVisc_Ns.clear();
m_coeffVisc_Ns.resize(m_nsp);
m_viscTempDepType_Ns.resize(m_nsp);
//The following two are not yet filled in LiquidTransportParams //for each species, assign viscosity model and coefficients
m_visc_Eij = tr.visc_Eij ; for (k = 0; k < m_nsp; k++) {
m_visc_Sij = tr.visc_Sij ; Cantera::LiquidTransportData &ltd = tr.LTData[k];
//specify temperature dependence
m_viscTempDepType_Ns[k] = ltd.model_viscosity;
//vector kentry corresponds to the k-th entry of m_coeffVisc_Ns
vector_fp &kentry = m_coeffVisc_Ns[k];
//save logarithm of pre-exponential for easier computation if ( m_viscTempDepType_Ns[k] == LTR_MODEL_CONSTANT
m_visc_logA.resize(m_nsp); || m_viscTempDepType_Ns[k] == LTR_MODEL_POLY ) {
for ( int i = 0; i < m_nsp; i++ ) kentry = ltd.viscCoeffs;
m_visc_logA[i] = log( m_visc_A[i] );
m_thermCond_A = tr.thermCond_A ; } else if ( m_viscTempDepType_Ns[k] == LTR_MODEL_ARRHENIUS ) {
m_thermCond_n = tr.thermCond_n ; kentry = ltd.viscCoeffs;
m_thermCond_Tact = tr.thermCond_Tact ; //for Arrhenius form, also carry the logarithm of the pre-exponential
kentry[3] = log( kentry[0] );
m_hydrodynamic_radius = tr.hydroRadius ; } else if ( m_viscTempDepType_Ns[k] == LTR_MODEL_NOTSET ) {
//we might be OK with viscosity not being set so
// this error is repeated in updateViscosity_T()
// and can be deleted from here if appropriate
throw CanteraError("LiquidTransport::initLiquid",
"Viscosity Model is not set for species "
+ m_thermo->speciesName(k)
+ " in the input file");
} else {
throw CanteraError("LiquidTransport::initLiquid",
"Viscosity Model for species "
+ m_thermo->speciesName(k)
+ " is not handled by this object");
}
}
/*
* Get the input Thermal Conductivities
*/
m_lambdaSpecies.resize(m_nsp);
m_coeffLambda_Ns.clear();
m_coeffLambda_Ns.resize(m_nsp);
m_lambdaTempDepType_Ns.resize(m_nsp);
//m_diffcoeffs = tr.diffcoeffs; //for each species, assign viscosity model and coefficients
for (k = 0; k < m_nsp; k++) {
Cantera::LiquidTransportData &ltd = tr.LTData[k];
//specify temperature dependence
m_lambdaTempDepType_Ns[k] = ltd.model_thermalCond;
//vector kentry corresponds to the k-th entry of m_coeffLambda_Ns
vector_fp &kentry = m_coeffLambda_Ns[k];
if ( m_lambdaTempDepType_Ns[k] == LTR_MODEL_CONSTANT
|| m_lambdaTempDepType_Ns[k] == LTR_MODEL_POLY ) {
kentry = ltd.thermalCondCoeffs;
} else if ( m_lambdaTempDepType_Ns[k] == LTR_MODEL_ARRHENIUS ) {
kentry = ltd.thermalCondCoeffs;
//for Arrhenius form, also carry the logarithm of the pre-exponential
kentry[3] = log( kentry[0] );
} else if ( m_lambdaTempDepType_Ns[k] == LTR_MODEL_NOTSET ) {
throw CanteraError("LiquidTransport::initLiquid",
"Thermal conductivity model is not set for species "
+ m_thermo->speciesName(k)
+ " in the input file");
} else {
throw CanteraError("LiquidTransport::initLiquid",
"Thermal conductivity model for species "
+ m_thermo->speciesName(k)
+ " is not handled by this object");
}
}
/*
* Get the input Hydrodynamic Radii
*/
m_hydrodynamic_radius.resize(m_nsp);
m_coeffRadius_Ns.clear();
m_coeffRadius_Ns.resize(m_nsp);
m_radiusTempDepType_Ns.resize(m_nsp);
//for each species, assign viscosity model and coefficients
for (k = 0; k < m_nsp; k++) {
Cantera::LiquidTransportData &ltd = tr.LTData[k];
//specify temperature dependence
m_radiusTempDepType_Ns[k] = ltd.model_hydroradius;
//vector kentry corresponds to the k-th entry of m_coeffRadius_Ns
vector_fp &kentry = m_coeffRadius_Ns[k];
if ( m_radiusTempDepType_Ns[k] == LTR_MODEL_CONSTANT
|| m_radiusTempDepType_Ns[k] == LTR_MODEL_POLY ) {
kentry = ltd.hydroRadiusCoeffs;
} else if ( m_radiusTempDepType_Ns[k] == LTR_MODEL_ARRHENIUS ) {
kentry = ltd.hydroRadiusCoeffs;
//for Arrhenius form, also carry the logarithm of the pre-exponential
kentry[3] = log( kentry[0] );
} else if ( m_radiusTempDepType_Ns[k] == LTR_MODEL_NOTSET ) {
throw CanteraError("LiquidTransport::initLiquid",
"Hydrodynamic radius model is not set for species "
+ m_thermo->speciesName(k)
+ " in the input file");
} else {
throw CanteraError("LiquidTransport::initLiquid",
"Hydrodynamic radius model for species "
+ m_thermo->speciesName(k)
+ " is not handled by this object");
}
}
/*
* Get the input Species Diffusivities
* Note that species diffusivities are not what is needed.
* Rather the Stefan Boltzmann interaction parameters are
* needed for the current model. This section may, therefore,
* be extraneous.
*/
// m_viscSpecies.resize(m_nsp);
m_coeffDiff_Ns.clear();
m_coeffDiff_Ns.resize(m_nsp);
m_diffTempDepType_Ns.resize(m_nsp);
//for each species, assign viscosity model and coefficients
for (k = 0; k < m_nsp; k++) {
Cantera::LiquidTransportData &ltd = tr.LTData[k];
//specify temperature dependence
if ( ltd.model_speciesDiffusivity >= 0
|| ltd.speciesDiffusivityCoeffs.size() > 0 ) {
cout << "Warning: diffusion coefficient data for "
<< m_thermo->speciesName(k)
<< endl
<< "in the input file is not used for LiquidTransport model."
<< endl
<< "LiquidTransport model uses hydrodynamicRadius, viscosity "
<< endl
<< "and the Stokes-Einstein equation."
<< endl;
}
}
m_mode = tr.mode_; m_mode = tr.mode_;
m_viscSpecies.resize(m_nsp); m_viscSpecies.resize(m_nsp);
m_logViscSpecies.resize(m_nsp); m_logViscSpecies.resize(m_nsp);
m_condSpecies.resize(m_nsp); m_lambdaSpecies.resize(m_nsp);
m_bdiff.resize(m_nsp, m_nsp); m_bdiff.resize(m_nsp, m_nsp);
m_molefracs.resize(m_nsp); m_molefracs.resize(m_nsp);
@ -237,22 +407,23 @@ namespace Cantera {
*/ */
doublereal LiquidTransport::viscosity() { doublereal LiquidTransport::viscosity() {
update_temp(); update_T();
update_conc(); update_C();
if (m_visc_mix_ok) return m_viscmix; if (m_visc_mix_ok) return m_viscmix;
// update m_viscSpecies[] if necessary // update m_viscSpecies[] if necessary
if (!m_visc_temp_ok) { if (!m_visc_temp_ok) {
updateViscosity_temp(); updateViscosity_T();
} }
if (!m_visc_conc_ok) { if (!m_visc_conc_ok) {
updateViscosities_conc(); updateViscosities_C();
} }
/* We still need to implement interaction parameters */ /* We still need to implement interaction parameters */
/* This constant viscosity model has no input */ /* This constant viscosity model has no input */
if (viscosityModel_ == LVISC_CONSTANT) { if (viscosityModel_ == LVISC_CONSTANT) {
err("constant viscosity not implemented for LiquidTransport."); err("constant viscosity not implemented for LiquidTransport.");
@ -272,15 +443,17 @@ namespace Cantera {
* ( m_visc_Sij(i,j) + m_visc_Eij(i,j) / m_temp ); * ( m_visc_Sij(i,j) + m_visc_Eij(i,j) / m_temp );
m_viscmix = exp( interaction ); m_viscmix = exp( interaction );
} else {
err("Unknown viscosity model in LiquidTransport::viscosity().");
} }
return m_viscmix; return m_viscmix;
} }
void LiquidTransport::getSpeciesViscosities(doublereal* visc) { void LiquidTransport::getSpeciesViscosities(doublereal* visc) {
update_temp(); update_T();
if (!m_visc_temp_ok) { if (!m_visc_temp_ok) {
updateViscosity_temp(); updateViscosity_T();
} }
copy(m_viscSpecies.begin(), m_viscSpecies.end(), visc); copy(m_viscSpecies.begin(), m_viscSpecies.end(), visc);
} }
@ -292,19 +465,23 @@ namespace Cantera {
void LiquidTransport::getBinaryDiffCoeffs(int ld, doublereal* d) { void LiquidTransport::getBinaryDiffCoeffs(int ld, doublereal* d) {
int i,j; int i,j;
update_temp(); if ( ld != m_nsp )
throw CanteraError("LiquidTransport::getBinaryDiffCoeffs",
"First argument does not correspond to number of species in model.\nDiff Coeff matrix may be misdimensioned");
update_T();
// if necessary, evaluate the binary diffusion coefficents // if necessary, evaluate the binary diffusion coefficents
// from the polynomial fits // from the polynomial fits
if (!m_diff_temp_ok) updateDiff_temp(); if (!m_diff_temp_ok) updateDiff_T();
doublereal pres = m_thermo->pressure();
doublereal rp = 1.0/pres;
for (i = 0; i < m_nsp; i++) for (i = 0; i < m_nsp; i++)
for (j = 0; j < m_nsp; j++) { for (j = 0; j < m_nsp; j++) {
d[ld*j + i] = rp * m_bdiff(i,j); d[ld*j + i] = m_bdiff(i,j);
} }
} }
//================================================================================================ //================================================================================================
// Get the electrical Mobilities (m^2/V/s). // Get the electrical Mobilities (m^2/V/s).
/* /*
@ -379,29 +556,32 @@ namespace Cantera {
} }
//================================================================================================ //================================================================================================
/****************** thermal conductivity **********************/ /****************** thermal conductivity **********************/
/* /*
* The thermal conductivity is computed from the following mixture rule: * The thermal conductivity is computed from the following mixture rule:
* \[ * \[
* \lambda = 0.5 \left( \sum_k X_k \lambda_k * \lambda = \left( \sum_k Y_k \lambda_k \right)
* + \frac{1}{\sum_k X_k/\lambda_k}\right) * \]
* \]
*/ */
doublereal LiquidTransport::thermalConductivity() { doublereal LiquidTransport::thermalConductivity() {
update_temp(); update_T();
update_conc(); update_C();
if (!m_cond_temp_ok) { if (!m_cond_temp_ok) {
updateCond_temp(); updateCond_T();
} }
if (!m_cond_mix_ok) { if (!m_cond_mix_ok) {
doublereal sum1 = 0.0, sum2 = 0.0;
for (int k = 0; k < m_nsp; k++) { // mass-fraction weighted thermal conductivity
sum1 += m_molefracs[k] * m_condSpecies[k]; {
sum2 += m_molefracs[k] / m_condSpecies[k]; doublereal sum1 = 0.0, sum2 = 0.0;
for (int k = 0; k < m_nsp; k++) {
sum1 += m_molefracs[k] * m_mw[k] * m_lambdaSpecies[k];
sum2 += m_molefracs[k] * m_mw[k] ;
}
m_lambda = sum1 / sum2 ;
} }
m_lambda = 0.5*(sum1 + 1.0/sum2);
m_cond_mix_ok = true; m_cond_mix_ok = true;
} }
@ -455,8 +635,8 @@ namespace Cantera {
void LiquidTransport::getSpeciesFluxesExt(int ldf, doublereal* fluxes) { void LiquidTransport::getSpeciesFluxesExt(int ldf, doublereal* fluxes) {
int n, k; int n, k;
update_temp(); update_T();
update_conc(); update_C();
getMixDiffCoeffs(DATA_PTR(m_spwork)); getMixDiffCoeffs(DATA_PTR(m_spwork));
@ -491,12 +671,12 @@ namespace Cantera {
*/ */
void LiquidTransport::getMixDiffCoeffs(doublereal* const d) { void LiquidTransport::getMixDiffCoeffs(doublereal* const d) {
update_temp(); update_T();
update_conc(); update_C();
// update the binary diffusion coefficients if necessary // update the binary diffusion coefficients if necessary
if (!m_diff_temp_ok) { if (!m_diff_temp_ok) {
updateDiff_temp(); updateDiff_T();
} }
int k, j; int k, j;
@ -534,20 +714,20 @@ namespace Cantera {
* This is called whenever a transport property is * This is called whenever a transport property is
* requested. * requested.
* The first task is to check whether the temperature has changed * The first task is to check whether the temperature has changed
* since the last call to update_temp(). * since the last call to update_T().
* If it hasn't then an immediate return is carried out. * If it hasn't then an immediate return is carried out.
* *
* @internal * @internal
*/ */
void LiquidTransport::update_temp() bool LiquidTransport::update_T()
{ {
// First make a decision about whether we need to recalculate // First make a decision about whether we need to recalculate
doublereal t = m_thermo->temperature(); doublereal t = m_thermo->temperature();
if (t == m_temp) return; if (t == m_temp) return false;
// Next do a reality check on temperature value // Next do a reality check on temperature value
if (t < 0.0) { if (t < 0.0) {
throw CanteraError("LiquidTransport::update_temp()", throw CanteraError("LiquidTransport::update_T()",
"negative temperature "+fp2str(t)); "negative temperature "+fp2str(t));
} }
@ -570,7 +750,7 @@ namespace Cantera {
m_visc_mix_ok = false; m_visc_mix_ok = false;
m_diff_mix_ok = false; m_diff_mix_ok = false;
// m_cond_mix_ok = false; (don't need it because a lower lvl flag is set // m_cond_mix_ok = false; (don't need it because a lower lvl flag is set
return true;
} }
@ -586,7 +766,7 @@ namespace Cantera {
* *
* @internal * @internal
*/ */
void LiquidTransport::update_conc() { bool LiquidTransport::update_C() {
// If the pressure has changed then the concentrations // If the pressure has changed then the concentrations
// have changed. // have changed.
doublereal pres = m_thermo->pressure(); doublereal pres = m_thermo->pressure();
@ -613,7 +793,7 @@ namespace Cantera {
concTot_tran_ *= concTot_; concTot_tran_ *= concTot_;
} }
if (qReturn) { if (qReturn) {
return; return false;
} }
// signal that concentration-dependent quantities will need to // signal that concentration-dependent quantities will need to
@ -625,6 +805,8 @@ namespace Cantera {
m_visc_mix_ok = false; m_visc_mix_ok = false;
m_diff_mix_ok = false; m_diff_mix_ok = false;
m_cond_mix_ok = false; m_cond_mix_ok = false;
return true;
} }
@ -680,71 +862,84 @@ namespace Cantera {
/************************************************************************* /*************************************************************************
* *
* methods to update temperature-dependent properties * methods to update species temperature-dependent properties
* *
*************************************************************************/ *************************************************************************/
/** /**
* Update the temperature-dependent parts of the mixture-averaged * Update the temperature-dependent parts of the species
* thermal conductivity. * thermal conductivity.
*/ */
void LiquidTransport::updateCond_temp() { void LiquidTransport::updateCond_T() {
int k;
/* for (k = 0; k < m_nsp; k++) {
if (m_mode == CK_Mode) { vector_fp &coeffk = m_coeffLambda_Ns[k];
for (k = 0; k < m_nsp; k++) {
m_condSpecies[k] = exp(m_condcoeffs[k]); if ( m_lambdaTempDepType_Ns[k] == LTR_MODEL_CONSTANT ) {
} m_lambdaSpecies[k] = coeffk[0] ;
} else {
for (k = 0; k < m_nsp; k++) { } else if ( m_lambdaTempDepType_Ns[k] == LTR_MODEL_ARRHENIUS ) {
m_condSpecies[k] = m_sqrt_t * m_condcoeffs[k]; //m_coeffLambda_Ns[k][0] holds A
//m_coeffLambda_Ns[k][1] holds n
//m_coeffLambda_Ns[k][2] holds Tact
//m_coeffLambda_Ns[k][3] holds log(A)
m_lambdaSpecies[k] = coeffk[0] * exp( coeffk[1] * m_logt
- coeffk[2] / m_temp );
} else if ( m_lambdaTempDepType_Ns[k] == LTR_MODEL_POLY ) {
m_lambdaSpecies[k] = coeffk[0]
+ coeffk[1] * m_temp
+ coeffk[2] * m_temp * m_temp
+ coeffk[3] * m_temp * m_temp * m_temp
+ coeffk[4] * m_temp * m_temp * m_temp * m_temp;
} else if ( m_lambdaTempDepType_Ns[k] == LTR_MODEL_NOTSET ) {
throw CanteraError("LiquidTransport::updateCond_T",
"Conductivity Model is not set for species "
+ m_thermo->speciesName(k)
+ " in the input file");
} else {
throw CanteraError("LiquidTransport::updateCond_T",
"Conductivity Model for species "
+ m_thermo->speciesName(k)
+ " is not handled by this object");
} }
} }
m_cond_temp_ok = true; m_cond_temp_ok = true;
m_cond_mix_ok = false; m_cond_mix_ok = false;
*/
} }
//! Update the StefanMaxwell interaction parameters.
/** /**
* Update the binary diffusion coefficients. These are evaluated * These are evaluated using the Stokes-Einstein
* from the polynomial fits at unit pressure (1 Pa). * relation from the viscosity and hydrodynamic radius.
*/ */
void LiquidTransport::updateDiff_temp() { void LiquidTransport::updateDiff_T() {
// evaluate binary diffusion coefficients at unit pressure double *viscSpec = new double(m_nsp);
double *radiusSpec = new double(m_nsp);
getSpeciesViscosities( viscSpec );
getSpeciesHydrodynamicRadius( radiusSpec );
/* int i,j;
if (m_mode == CK_Mode) { for (i = 0; i < m_nsp; i++)
for (i = 0; i < m_nsp; i++) { for (j = 0; j < m_nsp; j++) {
for (j = i; j < m_nsp; j++) { m_DiffCoeff_StefMax(i,j) = m_bdiff(i,j) = GasConstant * m_temp
m_bdiff(i,j) = exp(m_diffcoeffs[ic]); / ( 6.0 * Pi * radiusSpec[i] * viscSpec[j] ) ;
m_bdiff(j,i) = m_bdiff(i,j); cout << " D_ij = " << m_bdiff(i,j) << " for "
ic++; << m_thermo->speciesName(i) << ", "
} << m_thermo->speciesName(j) << endl;
} }
}
else {
for (i = 0; i < m_nsp; i++) {
for (j = i; j < m_nsp; j++) {
m_bdiff(i,j) = m_temp * m_sqrt_t*m_diffcoeffs[ic];
m_bdiff(j,i) = m_bdiff(i,j);
ic++;
}
}
}
m_diff_temp_ok = true; m_diff_temp_ok = true;
m_diff_mix_ok = false; m_diff_mix_ok = false;
*/
} }
/** //! Update the pure-species viscosities functional dependence on concentration.
* Update the pure-species viscosities. void LiquidTransport::updateViscosities_C() {
*/
void LiquidTransport::updateViscosities_conc() {
m_visc_conc_ok = true; m_visc_conc_ok = true;
} }
@ -755,20 +950,47 @@ namespace Cantera {
* weighting functions in the viscosity mixture rule. * weighting functions in the viscosity mixture rule.
* The flag m_visc_ok is set to true. * The flag m_visc_ok is set to true.
*/ */
void LiquidTransport::updateViscosity_temp() { void LiquidTransport::updateViscosity_T() {
int k; int k;
for (k = 0; k < m_nsp; k++) { for (k = 0; k < m_nsp; k++) {
m_logViscSpecies[k] = m_visc_logA[k] + m_visc_n[k] * m_logt vector_fp &coeffk = m_coeffVisc_Ns[k];
+ m_visc_Tact[k] / m_temp ;
m_viscSpecies[k] = exp( m_logViscSpecies[k] ); if ( m_viscTempDepType_Ns[k] == LTR_MODEL_CONSTANT ) {
m_logViscSpecies[k] = log( coeffk[0] );
m_viscSpecies[k] = coeffk[0] ;
} else if ( m_viscTempDepType_Ns[k] == LTR_MODEL_ARRHENIUS ) {
//m_coeffVisc_Ns[k][0] holds A
//m_coeffVisc_Ns[k][1] holds n
//m_coeffVisc_Ns[k][2] holds Tact
//m_coeffVisc_Ns[k][3] holds log(A)
m_logViscSpecies[k] = coeffk[3] + coeffk[1] * m_logt
- coeffk[2] / m_temp ;
m_viscSpecies[k] = exp( m_logViscSpecies[k] );
} else if ( m_viscTempDepType_Ns[k] == LTR_MODEL_POLY ) {
m_viscSpecies[k] = coeffk[0]
+ coeffk[1] * m_temp
+ coeffk[2] * m_temp * m_temp
+ coeffk[3] * m_temp * m_temp * m_temp
+ coeffk[4] * m_temp * m_temp * m_temp * m_temp;
m_logViscSpecies[k] = log( m_viscSpecies[k] );
} else if ( m_viscTempDepType_Ns[k] == LTR_MODEL_NOTSET ) {
throw CanteraError("LiquidTransport::updateViscosity_T",
"Viscosity Model is not set for species "
+ m_thermo->speciesName(k)
+ " in the input file");
} else {
throw CanteraError("LiquidTransport::updateViscosity_T",
"Viscosity Model for species "
+ m_thermo->speciesName(k)
+ " is not handled by this object");
}
m_visc_temp_ok = true;
m_visc_mix_ok = false;
} }
//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;
} }
@ -787,9 +1009,10 @@ namespace Cantera {
/* /*
* Update the concentrations in the mixture. * Update the concentrations and diffusion coefficients in the mixture.
*/ */
update_conc(); update_C();
if ( !m_diff_temp_ok ) updateDiff_T();
double T = m_thermo->temperature(); double T = m_thermo->temperature();

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

@ -34,9 +34,9 @@ namespace Cantera {
const int LVISC_INTERACTION = 1; const int LVISC_INTERACTION = 1;
const int LVISC_AVG_ENERGIES = 2; const int LVISC_AVG_ENERGIES = 2;
const int LDIFF_MIXDIFF_UNCORRECTED = 0; const int LDIFF_CONSTANT = 0;
const int LDIFF_MIXDIFF_FLUXCORRECTED = 1; const int LDIFF_ARHENNIUS = 1;
const int LDIFF_MULTICOMP_STEFANMAXWELL = 2; const int LDIFF_STOKES_EINSTEIN = 2;
@ -136,6 +136,9 @@ namespace Cantera {
class LiquidTransport : public Transport { class LiquidTransport : public Transport {
public: public:
typedef vector_fp Coeff_T_;
//! Default constructor. //! Default constructor.
/*! /*!
* This requires call to initLiquid(LiquidTransportParams& tr) * This requires call to initLiquid(LiquidTransportParams& tr)
@ -222,10 +225,18 @@ namespace Cantera {
*/ */
virtual void getSpeciesViscosities(doublereal* const visc); virtual void getSpeciesViscosities(doublereal* const visc);
//! Returns the hydrodynamic radius for all species
/*!
* The pure species viscosities are to be given in an Arrhenius
* form in accordance with activated-jump-process dominated transport.
*/
virtual void getSpeciesHydrodynamicRadius(doublereal* const radius);
//! Returns the binary diffusion coefficients //! Returns the binary diffusion coefficients
/*! /*!
* @param ld * @param ld number of species in system
* @param d * @param d vector of mixture diffusion coefficients
* units = m2 s-1. length = ld*ld = (number of species)^2
*/ */
virtual void getBinaryDiffCoeffs(const int ld, doublereal* const d); virtual void getBinaryDiffCoeffs(const int ld, doublereal* const d);
@ -237,14 +248,19 @@ namespace Cantera {
virtual void getMixDiffCoeffs(doublereal* const d); 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); virtual void getThermalDiffCoeffs(doublereal* const dt);
//! Return the thermal conductivity of the solution //! Return the thermal conductivity of the solution
/*! /*!
* The thermal conductivity is computed from the following mixture rule: * The thermal conductivity is computed from the following mixture rule:
* \f[ * \f[
* \lambda = 0.5 \left( \sum_k X_k \lambda_k * \lambda = \left( \sum_k Y_k \lambda_k \right)
* + \frac{1}{\sum_k X_k/\lambda_k}\right)
* \f] * \f]
* *
* Controlling update boolean = m_condmix_ok * Controlling update boolean = m_condmix_ok
@ -309,14 +325,84 @@ namespace Cantera {
*/ */
virtual void set_Grad_X(const doublereal* const grad_X); virtual void set_Grad_X(const doublereal* const grad_X);
virtual void update_Grad_lnAC();
//! Updates the internal value of the gradient of the logarithm of the
//! activity coefficients, which is used in the gradient of the chemical potential.
/*! The gradient of the chemical potential can be written in terms of
* gradient of the logarithm of the mole fraction times a correction
* associated with the gradient of the activity coefficient relative to
* that of the mole fraction. Specifically, the gradients of the
* logarithms of each are involved according to the formula
* \f[
* \nabla \mu_k = RT \nabla ( \ln X_k )
* \[ 1 + \nabla ( \ln \gamma_k ) / \nabla ( \ln X_k ) \]
* \f]
*
* The quantity within the square brackets is computed within
* this method.
*/
virtual void update_Grad_lnAC();
/**
* @param ndim The number of spatial dimensions (1, 2, or 3).
* @param grad_T The temperature gradient (ignored in this model).
* (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 fluxes Output of the diffusive mass fluxes
* Flat vector with the m_nsp in the inner loop.
* length = ldx * ndim
*
*
* 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: protected:
//! Handles the effects of changes in the Temperature, internally //! Returns true if temperature has changed,
//! within the object. //! in which case flags are set to recompute transport properties.
/*! /*!
* This is called whenever a transport property is * This is called whenever a transport property is requested.
* requested.
* The first task is to check whether the temperature has changed * The first task is to check whether the temperature has changed
* since the last call to update_T(). * since the last call to update_T().
* If it hasn't then an immediate return is carried out. * If it hasn't then an immediate return is carried out.
@ -326,10 +412,13 @@ namespace Cantera {
* part of all of the interfaces. * part of all of the interfaces.
* *
* @internal * @internal
*
* @return Returns true if the temperature has changed, and false otherwise
*/ */
virtual void update_temp(); virtual bool update_T();
//! Handles the effects of changes in the mixture concentration //! Returns true if mixture composition has changed,
//! in which case flags are set to recompute transport properties.
/*! /*!
* This is called for every interface call to check whether * This is called for every interface call to check whether
* the concentrations have changed. Concentrations change * the concentrations have changed. Concentrations change
@ -340,43 +429,47 @@ namespace Cantera {
* part of all of the interfaces. * part of all of the interfaces.
* *
* @internal * @internal
*
* @return Returns true if the mixture composition has changed, and false otherwise.
*/ */
virtual void update_conc(); virtual bool update_C();
public:
/**
* @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);
/**
* @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 getSpeciesFluxesExt(int ldf, doublereal* fluxes);
//! Solve the stefan_maxell equations for the diffusive fluxes. //! Solve the stefan_maxell equations for the diffusive fluxes.
void stefan_maxwell_solve(); void stefan_maxwell_solve();
//! 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 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
* m_visc_conc_ok 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: private:
//! Number of species in the mixture //! Number of species in the mixture
int m_nsp; int m_nsp;
@ -392,11 +485,17 @@ namespace Cantera {
*/ */
vector_fp m_mw; vector_fp m_mw;
//! Pure species viscosities in Arrhenius temperature-dependent form. //! Viscosity temperature dependence type
vector_fp m_visc_A; /*!
vector_fp m_visc_logA; //logarithm of coefficient * Types of temperature dependencies:
vector_fp m_visc_n; * 0 - Independent of temperature (only one implemented so far)
vector_fp m_visc_Tact; * 1 - extended arrhenius form
* 2 - polynomial in temperature form
*/
vector<LiquidTR_Model> m_viscTempDepType_Ns;
//! Pure species viscosities in temperature-dependent form.
std::vector<Coeff_T_> m_coeffVisc_Ns;
//! Molecular interaction energies associated with viscosity //! Molecular interaction energies associated with viscosity
/** /**
@ -412,22 +511,76 @@ namespace Cantera {
*/ */
DenseMatrix m_visc_Sij; DenseMatrix m_visc_Sij;
//! Pure species thermal conductivities in Arrhenius temperature-dependent form.
vector_fp m_thermCond_A;
vector_fp m_thermCond_n; //! Thermal conductivity temperature dependence type
vector_fp m_thermCond_Tact; /*!
* Types of temperature dependencies:
* 0 - Independent of temperature (only one implemented so far)
* 1 - extended arrhenius form
* 2 - polynomial in temperature form
*/
vector<LiquidTR_Model> m_lambdaTempDepType_Ns;
//! Pure species thermal conductivities in temperature-dependent form.
std::vector<Coeff_T_> m_coeffLambda_Ns;
//! Diffusion coefficient temperature dependence type
/*!
* Types of temperature dependencies:
* 0 - Independent of temperature (only one implemented so far)
* 1 - extended arrhenius form
* 2 - polynomial in temperature form
*/
vector<LiquidTR_Model> m_diffTempDepType_Ns;
//! Pure species diffusvities in temperature-dependent form.
std::vector<Coeff_T_> m_coeffDiff_Ns;
vector<bool> useHydroRadius_;
//!Hydrodynamic radius temperature dependence type
/*!
* Types of temperature dependencies:
* 0 - Independent of temperature
* 1 - extended arrhenius form
* 2 - polynomial in temperature form
*/
vector<LiquidTR_Model> m_radiusTempDepType_Ns;
//! Pure hydrodynamic radius in temperature-dependent form.
std::vector<Coeff_T_> m_coeffRadius_Ns;
//! Species hydrodynamic radius //! Species hydrodynamic radius
vector_fp m_hydrodynamic_radius; vector_fp m_hydrodynamic_radius;
//! 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 m_compositionDepType;
//! Polynomial coefficients of the binary diffusion coefficients //! Polynomial coefficients of the binary diffusion coefficients
/*! /*!
* These express the temperature dependendence of the * These express the temperature dependendence of the
* binary diffusivities. An overall pressure dependence is then * binary diffusivities. An overall pressure dependence is then
* added. * added.
*/ */
/*
vector<vector_fp> m_diffcoeffs; vector<vector_fp> m_diffcoeffs;
*/
//! Internal value of the gradient of the mole fraction vector //! Internal value of the gradient of the mole fraction vector
/*! /*!
@ -532,7 +685,7 @@ namespace Cantera {
* *
* controlling update boolean -> m_cond_temp_ok * controlling update boolean -> m_cond_temp_ok
*/ */
vector_fp m_condSpecies; vector_fp m_lambdaSpecies;
//! State of the mole fraction vector. //! State of the mole fraction vector.
int m_iStateMF; int m_iStateMF;
@ -587,7 +740,10 @@ namespace Cantera {
*/ */
doublereal concTot_tran_; doublereal concTot_tran_;
//! Mean molecular mass
doublereal meanMolecularWeight_; doublereal meanMolecularWeight_;
//! Density
doublereal dens_; doublereal dens_;
//! Local copy of the charge of each species //! Local copy of the charge of each species
@ -651,39 +807,13 @@ namespace Cantera {
//! Saved value of the mixture viscosity //! Saved value of the mixture viscosity
doublereal m_viscmix; doublereal m_viscmix;
// work space //! work space
/*!
* Length is equal to m_nsp
*/
vector_fp m_spwork; vector_fp m_spwork;
//! Internal Function
protected:
//! 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 updateViscosity_temp();
//! Update the temperature-dependent parts of the mixture-averaged
//! thermal conductivity.
void updateCond_temp();
//! Update the concentration parts of the viscosities
/*!
* Internal routine is run whenever the update_boolean
* m_visc_conc_ok is false. This routine will calculate
* internal values for the species viscosities.
*
* @internal
*/
void updateViscosities_conc();
//! Update the binary diffusion coefficients wrt T.
/*!
* These are evaluated
* from the polynomial fits at unit pressure (1 Pa).
*/
void updateDiff_temp();
private: private:
//! Boolean indicating that the top-level mixture viscosity is current //! Boolean indicating that the top-level mixture viscosity is current
@ -696,7 +826,7 @@ namespace Cantera {
bool m_visc_temp_ok; bool m_visc_temp_ok;
//! Flag to indicate that the pure species viscosities //! Flag to indicate that the pure species viscosities
//! are current wrt the temperature //! are current wrt the concentration
bool m_visc_conc_ok; bool m_visc_conc_ok;
//! Boolean indicating that mixture diffusion coeffs are current //! Boolean indicating that mixture diffusion coeffs are current
@ -712,18 +842,9 @@ namespace Cantera {
//! Boolean indicating that mixture conductivity is current //! Boolean indicating that mixture conductivity is current
bool m_cond_mix_ok; bool m_cond_mix_ok;
//! Mode for fitting the species viscosities //! Mode indicator for transport models -- currently unused.
/*!
* Either its CK_Mode or its cantera mode
* in CK_Mode visc is fitted to a polynomial
* in Cantera mode sqrt(visc) is fitted.
*/
int m_mode; int m_mode;
//! Internal storage for the diameter - diameter
//! species interactions
DenseMatrix m_diam;
//! Debugging flags //! Debugging flags
/*! /*!
* Turn on to get debugging information * Turn on to get debugging information
@ -736,8 +857,6 @@ namespace Cantera {
*/ */
int m_nDim; 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. * This probably indicates something is not yet implemented.