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

Implemented stefan maxwell equations with input of conductivity, cation mobility ratio, and self diffusion coefficients.

Browse source
This commit is contained in:
Christopher Lueth 2010-03-01 16:16:15 +00:00
parent 0d31b2f7a7
commit 60afcc5d57
8 changed files with 1930 additions and 314 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

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

@ -34,23 +34,38 @@ namespace Cantera {
LiquidTransport::LiquidTransport(thermo_t* thermo, int ndim) :
Transport(thermo, ndim),
m_nsp(0),
m_nBinInt(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_ionCondmix(-1.0),
m_mobRatMix(-1.0),
m_selfDiffMix(-1.0),
m_visc_mix_ok(false),
m_visc_temp_ok(false),
m_visc_conc_ok(false),
m_ionCond_mix_ok(false),
m_ionCond_temp_ok(false),
m_ionCond_conc_ok(false),
m_mobRat_mix_ok(false),
m_mobRat_temp_ok(false),
m_mobRat_conc_ok(false),
m_selfDiff_mix_ok(false),
m_selfDiff_temp_ok(false),
m_selfDiff_conc_ok(false),
m_radi_mix_ok(false),
m_radi_temp_ok(false),
m_radi_conc_ok(false),
m_diff_mix_ok(false),
m_diff_temp_ok(false),
m_cond_temp_ok(false),
m_cond_mix_ok(false),
m_lambda_temp_ok(false),
m_lambda_mix_ok(false),
m_mode(-1000),
m_debug(false),
m_nDim(1)
@ -61,23 +76,38 @@ namespace Cantera {
LiquidTransport::LiquidTransport(const LiquidTransport &right) :
Transport(),
m_nsp(0),
m_nBinInt(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_ionCondmix(-1.0),
m_mobRatMix(-1.0),
m_selfDiffMix(-1.0),
m_visc_mix_ok(false),
m_visc_temp_ok(false),
m_visc_conc_ok(false),
m_ionCond_mix_ok(false),
m_ionCond_temp_ok(false),
m_ionCond_conc_ok(false),
m_mobRat_mix_ok(false),
m_mobRat_temp_ok(false),
m_mobRat_conc_ok(false),
m_selfDiff_mix_ok(false),
m_selfDiff_temp_ok(false),
m_selfDiff_conc_ok(false),
m_radi_mix_ok(false),
m_radi_temp_ok(false),
m_radi_conc_ok(false),
m_diff_mix_ok(false),
m_diff_temp_ok(false),
m_cond_temp_ok(false),
m_cond_mix_ok(false),
m_lambda_temp_ok(false),
m_lambda_mix_ok(false),
m_mode(-1000),
m_debug(false),
m_nDim(1)
@ -95,10 +125,16 @@ namespace Cantera {
}
Transport::operator=(right);
m_nsp = right.m_nsp;
m_nBinInt = right.m_nBinInt;
m_tmin = right.m_tmin;
m_tmax = right.m_tmax;
m_mw = right.m_mw;
m_viscTempDep_Ns = right.m_viscTempDep_Ns;
m_ionCondTempDep_Ns = right.m_ionCondTempDep_Ns;
m_mobRatTempDep_Ns = right.m_mobRatTempDep_Ns;
m_mobRatTempDepIndex = right.m_mobRatTempDepIndex;
m_selfDiffTempDep_Ns = right.m_selfDiffTempDep_Ns;
m_selfDiffTempDepIndex = right.m_selfDiffTempDepIndex;
m_lambdaTempDep_Ns = right.m_lambdaTempDep_Ns;
m_diffTempDep_Ns = right.m_diffTempDep_Ns;
m_radiusTempDep_Ns = right.m_radiusTempDep_Ns;
@ -109,9 +145,19 @@ namespace Cantera {
m_Grad_mu = right.m_Grad_mu;
m_bdiff = right.m_bdiff;
m_viscSpecies = right.m_viscSpecies;
m_ionCondSpecies = right.m_ionCondSpecies;
m_mobRatSpecies = right.m_mobRatSpecies;
m_mobRatSpeciesIndex = right.m_mobRatSpeciesIndex;
m_selfDiffSpecies = right.m_selfDiffSpecies;
m_selfDiffSpeciesIndex = right.m_selfDiffSpeciesIndex;
m_hydrodynamic_radius = right.m_hydrodynamic_radius;
m_lambdaSpecies = right.m_lambdaSpecies;
m_viscMixModel = right.m_viscMixModel;
m_ionCondMixModel = right.m_ionCondMixModel;
m_mobRatMixModel = right.m_mobRatMixModel;
m_mobRatMixModelIndex = right.m_mobRatMixModelIndex;
m_selfDiffMixModel = right.m_selfDiffMixModel;
m_selfDiffMixModelIndex = right.m_selfDiffMixModelIndex;
m_lambdaMixModel = right.m_lambdaMixModel;
m_diffMixModel = right.m_diffMixModel;
m_iStateMF = -1;
@ -123,7 +169,6 @@ namespace Cantera {
m_actCoeff = right.m_actCoeff;
m_Grad_lnAC = right.m_Grad_lnAC;
m_chargeSpecies = right.m_chargeSpecies;
m_volume_spec = right.m_volume_spec;
m_B = right.m_B;
m_A = right.m_A;
m_temp = right.m_temp;
@ -132,17 +177,31 @@ namespace Cantera {
m_Vdiff = right.m_Vdiff;
m_lambda = right.m_lambda;
m_viscmix = right.m_viscmix;
m_ionCondmix = right.m_ionCondmix;
m_mobRatMix = right.m_mobRatMix;
m_mobRatMixIndex = right.m_mobRatMixIndex;
m_selfDiffMix = right.m_selfDiffMix;
m_selfDiffMixIndex = right.m_selfDiffMixIndex;
m_spwork = right.m_spwork;
m_visc_mix_ok = false;
m_visc_temp_ok = false;
m_visc_conc_ok = false;
m_ionCond_mix_ok = false;
m_ionCond_temp_ok = false;
m_ionCond_conc_ok = false;
m_mobRat_mix_ok = false;
m_mobRat_temp_ok = false;
m_mobRat_conc_ok = false;
m_selfDiff_mix_ok = false;
m_selfDiff_temp_ok = false;
m_selfDiff_conc_ok = false;
m_radi_mix_ok = false;
m_radi_temp_ok = false;
m_radi_conc_ok = false;
m_diff_mix_ok = false;
m_diff_temp_ok = false;
m_cond_temp_ok = false;
m_cond_mix_ok = false;
m_lambda_temp_ok = false;
m_lambda_mix_ok = false;
m_mode = right.m_mode;
m_debug = right.m_debug;
m_nDim = right.m_nDim;
@ -158,18 +217,33 @@ namespace Cantera {
LiquidTransport::~LiquidTransport() {
//These are constructed in TransportFactory::newLTP
for ( int k = 0; k < m_nsp; k++) {
if ( m_viscTempDep_Ns[k] ) delete m_viscTempDep_Ns[k];
if ( m_lambdaTempDep_Ns[k] ) delete m_lambdaTempDep_Ns[k];
if ( m_radiusTempDep_Ns[k] ) delete m_radiusTempDep_Ns[k];
if ( m_diffTempDep_Ns[k] ) delete m_diffTempDep_Ns[k];
}
//These are constructed in TransportFactory::newLTI
if ( m_viscMixModel ) delete m_viscMixModel;
if ( m_lambdaMixModel ) delete m_lambdaMixModel;
if ( m_diffMixModel ) delete m_diffMixModel;
//if ( m_radiusMixModel ) delete m_radiusMixModel;
//These are constructed in TransportFactory::newLTP
for ( int k = 0; k < m_nsp; k++) {
if ( m_viscTempDep_Ns[k] ) delete m_viscTempDep_Ns[k];
if ( m_ionCondTempDep_Ns[k] ) delete m_ionCondTempDep_Ns[k];
for ( int l=0;l < m_nsp; l++ ){
if ( m_selfDiffTempDep_Ns[l][k] ) delete m_selfDiffTempDep_Ns[l][k];
}
for ( int l=0;l < m_nBinInt; l++ ){
if ( m_mobRatTempDep_Ns[l][k] ) delete m_mobRatTempDep_Ns[l][k];
}
if ( m_lambdaTempDep_Ns[k] ) delete m_lambdaTempDep_Ns[k];
if ( m_radiusTempDep_Ns[k] ) delete m_radiusTempDep_Ns[k];
if ( m_diffTempDep_Ns[k] ) delete m_diffTempDep_Ns[k];
//These are constructed in TransportFactory::newLTI
if ( m_selfDiffMixModel[k] ) delete m_selfDiffMixModel[k];
}
for ( int k = 0; k < m_nBinInt; k++) {
if ( m_mobRatMixModel[k] ) delete m_mobRatMixModel[k];
}
if ( m_viscMixModel ) delete m_viscMixModel;
if ( m_ionCondMixModel ) delete m_ionCondMixModel;
if ( m_lambdaMixModel ) delete m_lambdaMixModel;
if ( m_diffMixModel ) delete m_diffMixModel;
//if ( m_radiusMixModel ) delete m_radiusMixModel;
}
@ -191,6 +265,7 @@ namespace Cantera {
m_thermo = tr.thermo;
m_velocityBasis = tr.velocityBasis_;
m_nsp = m_thermo->nSpecies();
m_nBinInt = m_nsp*(m_nsp-1)/2;
m_tmin = m_thermo->minTemp();
m_tmax = m_thermo->maxTemp();
@ -200,38 +275,77 @@ namespace Cantera {
m_thermo->molecularWeights().end(), m_mw.begin());
/*
* Get the input Viscosities
* Get the input Viscosities, and stuff
*/
m_viscSpecies.resize(m_nsp, 0.0);
m_viscTempDep_Ns.resize(m_nsp, 0);
m_ionCondSpecies.resize(m_nsp, 0.0);
m_ionCondTempDep_Ns.resize(m_nsp, 0);
m_mobRatTempDepIndex.resize(m_nBinInt);
m_mobRatTempDep_Ns.resize(m_nBinInt);
m_mobRatMixModel.resize(m_nBinInt);
m_mobRatMixModelIndex.resize(m_nBinInt);
m_mobRatSpeciesIndex.resize(m_nBinInt);
m_mobRatSpecies.resize(m_nBinInt, m_nsp, 0.0);
m_mobRatMix.resize(m_nBinInt,0.0);
m_mobRatMixIndex.resize(m_nBinInt);
m_selfDiffTempDepIndex.resize(m_nsp);
m_selfDiffTempDep_Ns.resize(m_nsp);
m_selfDiffMixModel.resize(m_nsp);
m_selfDiffMixModelIndex.resize(m_nsp);
m_selfDiffSpeciesIndex.resize(m_nsp);
m_selfDiffSpecies.resize(m_nsp, m_nsp, 0.0);
m_selfDiffMix.resize(m_nsp,0.0);
m_selfDiffMixIndex.resize(m_nsp);
for (k=0; k < m_nsp; k++){
m_selfDiffTempDep_Ns[k].resize(m_nsp, 0);
}
for (k=0; k < m_nBinInt; k++){
m_mobRatTempDep_Ns[k].resize(m_nsp, 0);
}
m_lambdaSpecies.resize(m_nsp, 0.0);
m_lambdaTempDep_Ns.resize(m_nsp, 0);
m_hydrodynamic_radius.resize(m_nsp, 0.0);
m_radiusTempDep_Ns.resize(m_nsp, 0);
//first populate mixing rules and indices
for (k = 0; k < m_nsp; k++) {
m_selfDiffMixModel[k] = tr.selfDiffusion[k];
m_selfDiffMixModelIndex[k] = tr.selfDiffIndex[k];
}
for (k = 0; k < m_nBinInt; k++) {
m_mobRatMixModel[k] = tr.mobilityRatio[k];
m_mobRatMixModelIndex[k] = tr.mobRatIndex[k];
}
//for each species, assign viscosity model and coefficients
for (k = 0; k < m_nsp; k++) {
Cantera::LiquidTransportData &ltd = tr.LTData[k];
m_viscTempDep_Ns[k] = ltd.viscosity;
}
/*
* Get the input Thermal Conductivities
*/
m_lambdaSpecies.resize(m_nsp, 0.0);
m_lambdaTempDep_Ns.resize(m_nsp, 0);
//for each species, assign thermal conductivity model
for (k = 0; k < m_nsp; k++) {
Cantera::LiquidTransportData &ltd = tr.LTData[k];
m_ionCondTempDep_Ns[k] = ltd.ionConductivity;
for (int j = 0; j < m_nBinInt; j++){
for (int l=0; l < m_nBinInt; l++){
if (m_mobRatMixModelIndex[j] == ltd.mobRatIndex[l]) {
m_mobRatTempDep_Ns[j][k] = ltd.mobilityRatio[l];
m_mobRatTempDepIndex[j] = ltd.mobRatIndex[l];
break;
}
}
}
for (int j = 0; j < (int) m_selfDiffMixModelIndex.size(); j++){
for (int l=0; l < (int) m_selfDiffMixModelIndex.size(); l++){
if (m_selfDiffMixModelIndex[j] == ltd.selfDiffIndex[l]) {
m_selfDiffTempDep_Ns[j][k] = ltd.selfDiffusion[l];
m_selfDiffTempDepIndex[j] = ltd.selfDiffIndex[l];
break;
}
}
}
m_lambdaTempDep_Ns[k] = ltd.thermalCond;
}
/*
* Get the input Hydrodynamic Radii
*/
m_hydrodynamic_radius.resize(m_nsp, 0.0);
m_radiusTempDep_Ns.resize(m_nsp, 0);
//for each species, assign model for hydrodynamic radius
for (k = 0; k < m_nsp; k++) {
Cantera::LiquidTransportData &ltd = tr.LTData[k];
m_radiusTempDep_Ns[k] = ltd.hydroRadius;
}
/*
* Get the input Species Diffusivities
* Note that species diffusivities are not what is needed.
@ -263,18 +377,28 @@ namespace Cantera {
* species diffusivity and hydrodynamics radius (perhaps not needed in the
* present class).
*/
m_viscMixModel = tr.viscosity;
tr.viscosity = 0;
m_ionCondMixModel = tr.ionConductivity;
tr.ionConductivity = 0;
//m_mobRatMixModel = tr.mobilityRatio;
m_lambdaMixModel = tr.thermalCond;
tr.thermalCond = 0;
m_diffMixModel = tr.speciesDiffusivity;
tr.speciesDiffusivity = 0;
m_bdiff.resize(m_nsp,m_nsp, 0.0);
//Don't really need to update this here.
//It is updated in updateDiff_T()
m_diffMixModel->getMatrixTransProp(m_bdiff);
m_diffMixModel->getMatrixTransProp( m_bdiff );
m_mode = tr.mode_;
@ -306,10 +430,19 @@ namespace Cantera {
m_visc_mix_ok = false;
m_visc_temp_ok = false;
m_visc_conc_ok = false;
m_ionCond_mix_ok = false;
m_ionCond_temp_ok = false;
m_ionCond_conc_ok = false;
m_mobRat_mix_ok = false;
m_mobRat_temp_ok = false;
m_mobRat_conc_ok = false;
m_selfDiff_mix_ok = false;
m_selfDiff_temp_ok = false;
m_selfDiff_conc_ok = false;
m_radi_temp_ok = false;
m_radi_conc_ok = false;
m_cond_temp_ok = false;
m_cond_mix_ok = false;
m_lambda_temp_ok = false;
m_lambda_mix_ok = false;
m_diff_temp_ok = false;
m_diff_mix_ok = false;
@ -357,6 +490,215 @@ namespace Cantera {
copy(m_viscSpecies.begin(), m_viscSpecies.end(), visc);
}
/****************** ionConductivity ******************************/
// Returns the ionic conductivity of the solution
/*
* The ionConductivity calculation is handled by subclasses of
* LiquidTranInteraction as specified in the input file.
* These in turn employ subclasses of LTPspecies to
* determine the individual species ionic conductivities.
*/
doublereal LiquidTransport:: ionConductivity() {
update_T();
update_C();
if (m_ionCond_mix_ok) return m_ionCondmix;
////// LiquidTranInteraction method
m_ionCondmix = m_ionCondMixModel->getMixTransProp( m_ionCondTempDep_Ns );
return m_ionCondmix;
/*
// update m_ionCondSpecies[] if necessary
if (!m_ionCond_temp_ok) {
updateIonConductivity_T();
}
if (!m_ionCond_conc_ok) {
updateIonConductivity_C();
}
*/
}
// Returns the pure species ionic conductivities for all species
/*
* The pure species ionic conductivities are evaluated using the
* appropriate subclasses of LTPspecies as specified in the
* input file.
*
* @param ionCond array of length "number of species"
* to hold returned ionic conductivities.
*/
void LiquidTransport::getSpeciesIonConductivity(doublereal* ionCond) {
update_T();
if (!m_ionCond_temp_ok) {
updateIonConductivity_T();
}
copy(m_ionCondSpecies.begin(), m_ionCondSpecies.end(), ionCond);
}
/****************** mobilityRatio ******************************/
// Returns the mobility ratios of the solution
/*
* The mobility ratio calculation is handled by subclasses of
* LiquidTranInteraction as specified in the input file.
* These in turn employ subclasses of LTPspecies to
* determine the individual species mobility ratios.
*/
void LiquidTransport:: mobilityRatio(vector_fp& mobRat, std::vector<std::string>& mobRatIndex) {
update_T();
update_C();
////// LiquidTranInteraction method
if (!m_mobRat_mix_ok){
for (int k = 0; k < m_nBinInt; k++){
if ( m_mobRatMixModelIndex[k] != m_mobRatTempDepIndex[k] )
throw CanteraError("LiquidTransport::mobilityRation","Mobility Ratio Indices Don't Match: Mixture vs. Species");
m_mobRatMix[k] = m_mobRatMixModel[k]->getMixTransProp( m_mobRatTempDep_Ns[k] );
m_mobRatMixIndex[k] = m_mobRatMixModelIndex[k];
}
}
for (int k = 0; k < m_nBinInt; k++){
mobRat[k] = m_mobRatMix[k];
mobRatIndex[k]= m_mobRatMixIndex[k];
}
}
void LiquidTransport:: mobilityRatio(doublereal* mobRat, std::vector<std::string>& mobRatIndex) {
update_T();
update_C();
////// LiquidTranInteraction method
if (!m_mobRat_mix_ok){
for (int k = 0; k < m_nBinInt; k++){
if ( m_mobRatMixModelIndex[k] != m_mobRatTempDepIndex[k] )
throw CanteraError("LiquidTransport::mobilityRatio","Mobility Ratio Indices Don't Match: Mixture vs. Species");
m_mobRatMix[k] = m_mobRatMixModel[k]->getMixTransProp( m_mobRatTempDep_Ns[k] );
m_mobRatMixIndex[k] = m_mobRatMixModelIndex[k];
}
}
for (int k = 0; k < m_nBinInt; k++){
mobRat[k] = m_mobRatMix[k];
mobRatIndex[k]= m_mobRatMixIndex[k];
}
}
// Returns the pure species mobility ratios for all species
/*
* The pure species mobility ratios are evaluated using the
* appropriate subclasses of LTPspecies as specified in the
* input file.
*
* @param mobRat array of length "number of species"
* to hold returned mobility ratio.
*/
void LiquidTransport::getSpeciesMobilityRatio(DenseMatrix& mobRat, std::vector<std::string>& mobRatIndex) {
update_T();
if (!m_mobRat_temp_ok) {
updateMobilityRatio_T();
}
mobRat = m_mobRatSpecies;
for (int k = 0; k < m_nBinInt; k++)
mobRatIndex[k] = m_mobRatSpeciesIndex[k];
}
void LiquidTransport::getSpeciesMobilityRatio(doublereal** mobRat, std::vector<std::string>& mobRatIndex) {
update_T();
if (!m_mobRat_temp_ok) {
updateMobilityRatio_T();
}
for (int k=0; k<m_nBinInt; k++){
for (int j=0; j < m_nsp; j++)
mobRat[k][j] = m_mobRatSpecies(k,j);
mobRatIndex[k] = m_mobRatSpeciesIndex[k];
}
}
/****************** SelfDiffusion ******************************/
// Returns the mobility ratios of the solution
/*
* The mobility ratio calculation is handled by subclasses of
* LiquidTranInteraction as specified in the input file.
* These in turn employ subclasses of LTPspecies to
* determine the individual species mobility ratios.
*/
void LiquidTransport:: selfDiffusion(vector_fp& selfDiff, std::vector<std::string>& selfDiffIndex) {
update_T();
update_C();
////// LiquidTranInteraction method
if (!m_selfDiff_mix_ok){
for (int k = 0; k < m_nsp; k++){
if ( m_selfDiffMixModelIndex[k] != m_selfDiffTempDepIndex[k] )
throw CanteraError("LiquidTransport::selfDiffusion","Self Diffusion Indices Don't Match: Mixture vs. Species");
m_selfDiffMix[k] = m_selfDiffMixModel[k]->getMixTransProp( m_selfDiffTempDep_Ns[k] );
m_selfDiffMixIndex[k] = m_selfDiffMixModelIndex[k];
}
}
for (int k = 0; k < m_nsp; k++){
selfDiff[k] = m_selfDiffMix[k];
selfDiffIndex[k]= m_selfDiffMixIndex[k];
}
}
void LiquidTransport:: selfDiffusion(double* selfDiff, std::vector<std::string>& selfDiffIndex) {
update_T();
update_C();
////// LiquidTranInteraction method
if (!m_selfDiff_mix_ok){
for (int k = 0; k < m_nsp; k++){
if ( m_selfDiffMixModelIndex[k] != m_selfDiffTempDepIndex[k] )
throw CanteraError("LiquidTransport::selfDiffusion","Self Diffusion Indices Don't Match: Mixture vs. Species");
m_selfDiffMix[k] = m_selfDiffMixModel[k]->getMixTransProp( m_selfDiffTempDep_Ns[k] );
m_selfDiffMixIndex[k] = m_selfDiffMixModelIndex[k];
}
}
for (int k = 0; k < m_nsp; k++){
selfDiff[k] = m_selfDiffMix[k];
selfDiffIndex[k]= m_selfDiffMixIndex[k];
}
}
// Returns the pure species self diffusion for all species
/*
* The pure species self diffusion coeffs are evaluated using the
* appropriate subclasses of LTPspecies as specified in the
* input file.
*
* @param selfDiff array of size "number of species"^2
* to hold returned self diffusion.
*/
void LiquidTransport::getSpeciesSelfDiffusion(DenseMatrix& selfDiff, std::vector<std::string>& selfDiffIndex) {
update_T();
if (!m_selfDiff_temp_ok) {
updateSelfDiffusion_T();
}
selfDiff = m_selfDiffSpecies;
for (int k = 0; k < m_nsp; k++)
selfDiffIndex[k] = m_selfDiffSpeciesIndex[k];
}
void LiquidTransport::getSpeciesSelfDiffusion(doublereal** selfDiff, std::vector<std::string>& selfDiffIndex) {
update_T();
if (!m_selfDiff_temp_ok) {
updateSelfDiffusion_T();
}
for (int k=0; k<m_nsp; k++){
for (int j=0; j < m_nsp; j++)
selfDiff[k][j] = m_selfDiffSpecies(k,j);
selfDiffIndex[k] = m_selfDiffSpeciesIndex[k];
}
}
//===============================================================
// Returns the hydrodynamic radius for all species
/*
@ -390,7 +732,7 @@ namespace Cantera {
update_T();
update_C();
if (!m_cond_mix_ok) {
if (!m_lambda_mix_ok) {
m_lambda = m_lambdaMixModel->getMixTransProp( m_lambdaTempDep_Ns );
m_cond_mix_ok = true;
}
@ -846,11 +1188,6 @@ namespace Cantera {
void LiquidTransport::getSpeciesVdiffExt(int ldf, doublereal* Vdiff) {
int n, k;
update_T();
update_C();
update_Grad_lnAC();
stefan_maxwell_solve();
for (n = 0; n < m_nDim; n++) {
@ -878,11 +1215,6 @@ namespace Cantera {
void LiquidTransport::getSpeciesFluxesExt(int ldf, doublereal* fluxes) {
int n, k;
update_T();
update_C();
update_Grad_lnAC();
stefan_maxwell_solve();
for (n = 0; n < m_nDim; n++) {
@ -926,11 +1258,6 @@ namespace Cantera {
*/
void LiquidTransport::getMixDiffCoeffs(doublereal* const d) {
update_T();
update_C();
update_Grad_lnAC();
stefan_maxwell_solve();
for ( int n = 0; n < m_nDim; n++) {
@ -975,19 +1302,28 @@ namespace Cantera {
// temperature has changed so temp flags are flipped
m_visc_temp_ok = false;
m_ionCond_temp_ok = false;
m_mobRat_temp_ok = false;
m_selfDiff_temp_ok = false;
m_radi_temp_ok = false;
m_diff_temp_ok = false;
m_lambda_temp_ok = false;
// temperature has changed, so polynomial temperature
// interpolations will need to be reevaluated.
// This means that many concentration
m_visc_conc_ok = false;
m_cond_temp_ok = false;
m_ionCond_conc_ok = false;
m_mobRat_conc_ok = false;
m_selfDiff_conc_ok = false;
// Mixture stuff needs to be evaluated
m_visc_mix_ok = false;
m_ionCond_mix_ok = false;
m_mobRat_mix_ok = false;
m_selfDiff_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_lambda_mix_ok = false; //(don't need it because a lower lvl flag is set
return true;
}
@ -1040,11 +1376,17 @@ namespace Cantera {
// be recomputed before use, and update the local mole
// fractions.
m_visc_conc_ok = false;
m_ionCond_conc_ok = false;
m_mobRat_conc_ok = false;
m_selfDiff_conc_ok = false;
// Mixture stuff needs to be evaluated
m_visc_mix_ok = false;
m_ionCond_mix_ok = false;
m_mobRat_mix_ok = false;
m_selfDiff_mix_ok = false;
m_diff_mix_ok = false;
m_cond_mix_ok = false;
m_lambda_mix_ok = false;
return true;
}
@ -1067,8 +1409,8 @@ namespace Cantera {
for (k = 0; k < m_nsp; k++) {
m_lambdaSpecies[k] = m_lambdaTempDep_Ns[k]->getSpeciesTransProp() ;
}
m_cond_temp_ok = true;
m_cond_mix_ok = false;
m_lambda_temp_ok = true;
m_lambda_mix_ok = false;
}
@ -1112,6 +1454,79 @@ namespace Cantera {
}
//! Update the pure-species ionic conductivities functional dependence on concentration.
void LiquidTransport::updateIonConductivity_C() {
m_ionCond_conc_ok = true;
}
/**
* Updates the array of pure species ionic conductivities internally
* using calls to the appropriate LTPspecies subclass.
* The flag m_ionCond_ok is set to true.
*/
void LiquidTransport::updateIonConductivity_T() {
int k;
for (k = 0; k < m_nsp; k++) {
m_ionCondSpecies[k] = m_ionCondTempDep_Ns[k]->getSpeciesTransProp() ;
}
m_ionCond_temp_ok = true;
m_ionCond_mix_ok = false;
}
//! Update the pure-species mobility ratios functional dependence on concentration.
void LiquidTransport::updateMobilityRatio_C() {
m_mobRat_conc_ok = true;
}
/**
* Updates the array of pure species mobility ratios internally
* using calls to the appropriate LTPspecies subclass.
* The flag m_mobRat_ok is set to true.
*/
void LiquidTransport::updateMobilityRatio_T() {
int k;
int j;
for (k = 0; k < m_nBinInt; k++) {
for (j = 0; j < m_nsp; j++) {
m_mobRatSpecies(k,j) = m_mobRatTempDep_Ns[k][j]->getSpeciesTransProp() ;
}
m_mobRatSpeciesIndex[k] = m_mobRatTempDepIndex[k];
}
m_mobRat_temp_ok = true;
m_mobRat_mix_ok = false;
}
//! Update the pure-species self diffusion functional dependence on concentration.
void LiquidTransport::updateSelfDiffusion_C() {
m_selfDiff_conc_ok = true;
}
/**
* Updates the array of pure species self diffusion internally
* using calls to the appropriate LTPspecies subclass.
* The flag m_selfDiff_ok is set to true.
*/
void LiquidTransport::updateSelfDiffusion_T() {
int k;
int j;
for (k = 0; k < m_nBinInt; k++) {
for (j = 0; j < m_nsp; j++) {
m_selfDiffSpecies(k,j) = m_selfDiffTempDep_Ns[k][j]->getSpeciesTransProp() ;
}
m_selfDiffSpeciesIndex[k] = m_selfDiffTempDepIndex[k];
}
m_selfDiff_temp_ok = true;
m_selfDiff_mix_ok = false;
}
//! Update the pure-species viscosities functional dependence on concentration.
void LiquidTransport::updateHydrodynamicRadius_C() {
m_radi_conc_ok = true;
@ -1164,14 +1579,23 @@ namespace Cantera {
int k;
vector_fp grad_lnAC(m_nsp);
m_thermo->getdlnActCoeffdlnC( DATA_PTR(grad_lnAC) );
doublereal grad_T;
vector_fp grad_lnAC(m_nsp), grad_X(m_nsp);
// IonsFromNeutralVPSSTP * tempIons = dynamic_cast<IonsFromNeutralVPSSTP *> m_thermo;
//MargulesVPSSTP * tempMarg = dynamic_cast<MargulesVPSSTP *> (tempIons->neutralMoleculePhase_);
for (k = 0; k < m_nsp; k++) {
m_Grad_lnAC[k] = grad_lnAC[k];
//m_thermo->getdlnActCoeffdlnX( DATA_PTR(grad_lnAC) );
for (k = 0; k < m_nDim; k++ ) {
grad_T = m_Grad_T[k];
grad_X.assign(m_Grad_X.begin()+m_nsp*k,m_Grad_X.begin()+m_nsp*(k+1));
m_thermo->getdlnActCoeff( grad_T, DATA_PTR(grad_X), DATA_PTR(grad_lnAC) );
for ( int i = 0; i < m_nsp; i++ )
grad_lnAC[i] += grad_X[i]/m_molefracs[i];
copy(grad_lnAC.begin(),grad_lnAC.end(),m_Grad_lnAC.begin()+m_nsp*k);
// std::cout << k << " m_Grad_lnAC = " << m_Grad_lnAC[k] << std::endl;
}
return;
}
@ -1220,10 +1644,13 @@ namespace Cantera {
//! grab a local copy of the molecular weights
const vector_fp& M = m_thermo->molecularWeights();
//! grad a local copy of the ion molar volume (inverse total ion concentration)
const doublereal vol = m_thermo->molarVolume();
/*
* Update the concentrations and diffusion coefficients in the mixture.
* Update the temperature, concentrations and diffusion coefficients in the mixture.
*/
update_T();
update_C();
if ( !m_diff_temp_ok ) updateDiff_T();
@ -1260,12 +1687,11 @@ namespace Cantera {
*
*/
for (i = 0; i < m_nsp; i++) {
double xi_denom = m_molefracs_tran[i];
for (a = 0; a < m_nDim; a++) {
m_Grad_mu[a*m_nsp + i] =
m_chargeSpecies[i] * Faraday * m_Grad_V[a]
//+ (m_volume_spec[i] - M[i]/dens_) * m_Grad_P[a]
+ GasConstant * T * m_Grad_X[a*m_nsp+i] * ( 1.0 + m_Grad_lnAC[i] ) / xi_denom;
+ GasConstant * T * m_Grad_lnAC[a*m_nsp+i];
}
}
@ -1286,8 +1712,11 @@ namespace Cantera {
* Just for Note, m_A(i,j) refers to the ith row and jth column.
* They are still fortran ordered, so that i varies fastest.
*/
double condSum1, condSum2;
switch (m_nDim) {
case 1: /* 1-D approximation */
m_B(0,0) = 0.0;
//equation for the reference velocity
for (j = 0; j < m_nsp; j++) {
@ -1298,9 +1727,10 @@ namespace Cantera {
else if ( ( m_velocityBasis >= 0 )
&& ( m_velocityBasis < m_nsp ) )
// use species number m_velocityBasis as reference velocity
if ( m_velocityBasis == j ) m_A(0,j) = 1.0;
else
throw CanteraError("LiquidTransport::stefan_maxwell_solve",
if ( m_velocityBasis == j ) m_A(0,j) = 1.0;
else m_A(0,j) = 0.0;
else
throw CanteraError("LiquidTransport::stefan_maxwell_solve",
"Unknown reference velocity provided.");
}
for (i = 1; i < m_nsp; i++){
@ -1320,7 +1750,30 @@ namespace Cantera {
//! invert and solve the system Ax = b. Answer is in m_B
solve(m_A, m_B);
/*
condSum2 = m_chargeSpecies[1]*m_chargeSpecies[1]*m_molefracs_tran[1]/m_bdiff(2,3) +
m_chargeSpecies[2]*m_chargeSpecies[2]*m_molefracs_tran[2]/m_bdiff(1,3) +
m_chargeSpecies[3]*m_chargeSpecies[3]*m_molefracs_tran[3]/m_bdiff(1,2);
condSum1 = m_molefracs_tran[1]/m_bdiff(1,2)/m_bdiff(1,3) +
m_molefracs_tran[2]/m_bdiff(2,3)/m_bdiff(1,2) +
m_molefracs_tran[3]/m_bdiff(1,3)/m_bdiff(2,3);
condSum2 = condSum2/condSum1*Faraday*Faraday/GasConstant/T/vol;
*/
condSum1 = 0;
for (i = 0; i < m_nsp; i++){
condSum1 -= Faraday*m_chargeSpecies[i]*m_B(i,0)*m_molefracs_tran[i]/vol;
}
/*
Check Mobility Ratio of Cations
cout << "mobility ratio = " << m_chargeSpecies[1]*(m_B(1,0)-m_B(2,0))/m_chargeSpecies[0]/(m_B(0,0)-m_B(2,0)) << endl;
*/
// cout << condSum1 << " = " << condSum2 << endl;
break;
case 2: /* 2-D approximation */
m_B(0,0) = 0.0;
@ -1335,8 +1788,9 @@ namespace Cantera {
&& ( m_velocityBasis < m_nsp ) )
// use species number m_velocityBasis as reference velocity
if ( m_velocityBasis == j ) m_A(0,j) = 1.0;
else
throw CanteraError("LiquidTransport::stefan_maxwell_solve",
else m_A(0,j) = 0.0;
else
throw CanteraError("LiquidTransport::stefan_maxwell_solve",
"Unknown reference velocity provided.");
}
for (i = 1; i < m_nsp; i++){
@ -1375,8 +1829,9 @@ namespace Cantera {
&& ( m_velocityBasis < m_nsp ) )
// use species number m_velocityBasis as reference velocity
if ( m_velocityBasis == j ) m_A(0,j) = 1.0;
else
throw CanteraError("LiquidTransport::stefan_maxwell_solve",
else m_A(0,j) = 0.0;
else
throw CanteraError("LiquidTransport::stefan_maxwell_solve",
"Unknown reference velocity provided.");
}
for (i = 1; i < m_nsp; i++){

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

@ -68,7 +68,7 @@ namespace Cantera {
* \f]
* for mass fraction \f$ Y_i \f$. For mole-averaged velocities
* \f[
* \sum_{i} X_i \vec{V_i} = 0
* \sum_{i} X_i \vec{V_i} = 0
* \f]
* for mole fraction \f$ X_i \f$.
*
@ -165,7 +165,7 @@ namespace Cantera {
* determine the individual species viscosities.
*/
virtual doublereal viscosity();
//! Returns the pure species viscosities for all species
/*!
* The pure species viscosities are evaluated using the
@ -177,6 +177,70 @@ namespace Cantera {
*/
virtual void getSpeciesViscosities(doublereal* const visc);
//! Returns the ionic conductivity of the solution
/*!
* The ionic conductivity calculation is handled by subclasses of
* LiquidTranInteraction as specified in the input file.
* These in turn employ subclasses of LTPspecies to
* determine the individual species ionic conductivities.
*/
virtual doublereal ionConductivity();
//! Returns the pure species ionic conductivities for all species
/*!
* The pure species ionic conductivities are evaluated using the
* appropriate subclasses of LTPspecies as specified in the
* input file.
*
* @param visc array of length "number of species"
* to hold returned ionic conductivities.
*/
virtual void getSpeciesIonConductivity(doublereal* const ionCond);
//! Returns the mobility ratio of the solution
/*!
* The mobility ratio calculation is handled by subclasses of
* LiquidTranInteraction as specified in the input file.
* These in turn employ subclasses of LTPspecies to
* determine the individual species mobility ratios.
*/
virtual void mobilityRatio(vector_fp& mobRat, std::vector<std::string>& mobRatIndex);
virtual void mobilityRatio(double* mobRat, std::vector<std::string>& mobRatIndex);
//! Returns the pure species mobility ratios for all species
/*!
* The pure species mobility ratios are evaluated using the
* appropriate subclasses of LTPspecies as specified in the
* input file.
*
* @param mobRat array of length "number of species"
* to hold returned mobility ratios.
*/
virtual void getSpeciesMobilityRatio(DenseMatrix& mobRat, std::vector<std::string>& mobRatIndex);
virtual void getSpeciesMobilityRatio(double** mobRat, std::vector<std::string>& mobRatIndex);
//! Returns the self diffusion coefficients in the solution
/*!
* The self diffusion calculation is handled by subclasses of
* LiquidTranInteraction as specified in the input file.
* These in turn employ subclasses of LTPspecies to
* determine the individual species self diffusion coeffs.
*/
virtual void selfDiffusion(vector_fp& selfDiff, std::vector<std::string>& selfDiffIndex);
virtual void selfDiffusion(double* selfDiff, std::vector<std::string>& selfDiffIndex);
//! Returns the pure species self diffusion in solution of each species
/*!
* The pure species molar volumes are evaluated using the
* appropriate subclasses of LTPspecies as specified in the
* input file.
*
* @param selfDiff array of length "number of species"
* to hold returned self diffusion coeffs.
*/
virtual void getSpeciesSelfDiffusion(DenseMatrix& selfDiff, std::vector<std::string>& selfDiffIndex);
virtual void getSpeciesSelfDiffusion(double** selfDiff, std::vector<std::string>& selfDiffIndex);
//! Returns the hydrodynamic radius for all species
/*!
* The species hydrodynamic radii are evaluated using the
@ -711,6 +775,25 @@ namespace Cantera {
*/
void updateViscosity_T();
//! Update the temperature-dependent ionic conductivity terms
//! for each species internally
/*!
* The flag m_ionCond_temp_ok is set to true.
*/
void updateIonConductivity_T();
//! Updates the array of pure species mobility ratios internally.
/*!
* The flag m_mobRat_ok is set to true.
*/
void updateMobilityRatio_T();
//! Updates the array of pure species self diffusion coeffs internally.
/*!
* The flag m_selfDiff_ok is set to true.
*/
void updateSelfDiffusion_T();
//! Update the temperature-dependent hydrodynamic radius terms
//! for each species internally
/*!
@ -731,7 +814,37 @@ namespace Cantera {
* @internal
*/
void updateViscosities_C();
//! Update the concentration parts of the ionic conductivity
/*!
* Internal routine is run whenever the update_boolean
* m_ionCond_conc_ok is false. Currently there is no concentration
* dependence for the pure species ionic conductivity.
*
* @internal
*/
void updateIonConductivity_C();
//! Update the concentration parts of the mobility ratio
/*!
* Internal routine is run whenever the update_boolean
* m_mobRat_conc_ok is false. Currently there is no concentration
* dependence for the pure species mobility ratio.
*
* @internal
*/
void updateMobilityRatio_C();
//! Update the concentration parts of the self diffusion
/*!
* Internal routine is run whenever the update_boolean
* m_selfDiff_conc_ok is false. Currently there is no concentration
* dependence for the pure species self diffusion.
*
* @internal
*/
void updateSelfDiffusion_C();
//! Update the concentration dependence of the hydrodynamic radius
/*!
* Internal routine is run whenever the update_boolean
@ -751,6 +864,7 @@ namespace Cantera {
//! Number of species in the mixture
int m_nsp;
int m_nBinInt;
//! Minimum temperature applicable to the transport property eval
doublereal m_tmin;
@ -782,6 +896,65 @@ namespace Cantera {
*/
LiquidTranInteraction *m_viscMixModel;
//! Ionic conductivity for each species expressed as an appropriate subclass
//! of LTPspecies
/*!
* These subclasses of LTPspecies evaluate the species-specific
* transport properties according to the parameters parsed in
* TransportFactory::getLiquidSpeciesTransportData().
*/
std::vector<LTPspecies*> m_ionCondTempDep_Ns;
//! Ionic Conductivity of the mixture expressed as a subclass of
//! LiquidTranInteraction
/*!
* These subclasses of LiquidTranInteraction evaluate the
* mixture transport properties according to the parameters parsed in
* TransportFactory::getLiquidInteractionsTransportData().
*/
LiquidTranInteraction *m_ionCondMixModel;
//! Mobility ratio for each species expressed as an appropriate subclass
//! of LTPspecies
/*!
* These subclasses of LTPspecies evaluate the species-specific
* transport properties according to the parameters parsed in
* TransportFactory::getLiquidSpeciesTransportData().
*/
typedef std::vector<LTPspecies*> LTPvector;
std::vector<LTPvector> m_mobRatTempDep_Ns;
std::vector<std::string> m_mobRatTempDepIndex;
//! Mobility ratio of the mixture expressed as a subclass of
//! LiquidTranInteraction
/*!
* These subclasses of LiquidTranInteraction evaluate the
* mixture transport properties according to the parameters parsed in
* TransportFactory::getLiquidInteractionsTransportData().
*/
std::vector<LiquidTranInteraction*> m_mobRatMixModel;
std::vector<std::string> m_mobRatMixModelIndex;
//! Self Diffusion for each species expressed as an appropriate subclass
//! of LTPspecies
/*!
* These subclasses of LTPspecies evaluate the species-specific
* transport properties according to the parameters parsed in
* TransportFactory::getLiquidSpeciesTransportData().
*/
std::vector<LTPvector> m_selfDiffTempDep_Ns;
std::vector<std::string> m_selfDiffTempDepIndex;
//! Self Diffusion of the mixture expressed as a subclass of
//! LiquidTranInteraction
/*!
* These subclasses of LiquidTranInteraction evaluate the
* mixture transport properties according to the parameters parsed in
* TransportFactory::getLiquidInteractionsTransportData().
*/
std::vector<LiquidTranInteraction*> m_selfDiffMixModel;
std::vector<std::string> m_selfDiffMixModelIndex;
//! Thermal conductivity for each species expressed as an
//! appropriate subclass of LTPspecies
/*!
@ -872,19 +1045,24 @@ namespace Cantera {
vector_fp m_Grad_X;
//! Gradient of the logarithm of the activity coefficients
//! with respect to the logarithm of the mole fraction, plus one.
/*!
* This quantity appears in the gradient of the chemical potential.
* It multiplies the gradient of the mole fraction, and in this way
* serves to "modify" the diffusion coefficient.
*
* m_Grad_lnAC[k] = \partial \left[ \ln ( \gamma_i ) \right]
* / \partial \left[ \ln ( \X_i ) \right]
* m_Grad_lnAC[k] = \nabla \ln ( \gamma_i ) + \nabla \ln ( \X_i )
*
* Note that where "mole fraction" is used here, whatever
* concentration-related variable applies, so that if
* molality is the concentration variable, the gradient of the
* activity coefficient should be with respect to the molality.
* 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_nDimm_
*
* m_Grad_X[n*m_nsp + k]
*
*/
vector_fp m_Grad_lnAC;
@ -964,6 +1142,41 @@ namespace Cantera {
*/
vector_fp m_viscSpecies;
//! Internal value of the species ionic conductivities
/*!
* Ionic conductivity of the species evaluated using subclass of LTPspecies
* held in m_ionCondTempDep_Ns.
*
* Length = number of species
*
* controlling update boolean -> m_ionCond_temp_ok
*/
vector_fp m_ionCondSpecies;
//! Internal value of the species mobility ratios
/*!
* Mobility ratio of the species evaluated using subclass of LTPspecies
* held in m_mobRatTempDep_Ns.
*
* Length = number of species
*
* controlling update boolean -> m_mobRat_temp_ok
*/
DenseMatrix m_mobRatSpecies;
std::vector<std::string> m_mobRatSpeciesIndex;
//! Internal value of the species self diffusion coefficients
/*!
* Self diffusion of the species evaluated using subclass of LTPspecies
* held in m_selfDiffTempDep_Ns.
*
* Length = number of species
*
* controlling update boolean -> m_selfDiff_temp_ok
*/
DenseMatrix m_selfDiffSpecies;
std::vector<std::string> m_selfDiffSpeciesIndex;
//! Internal value of the species individual thermal conductivities
/*!
* Thermal conductivities of the species evaluated using subclass
@ -1093,6 +1306,19 @@ namespace Cantera {
//! Saved value of the mixture viscosity
doublereal m_viscmix;
//! Saved value of the mixture ionic conductivity
doublereal m_ionCondmix;
//! Saved values of the mixture mobility ratios
vector_fp m_mobRatMix;
//! Saved species index of the mixture correlated to mobility ratios
std::vector<std::string> m_mobRatMixIndex;
//! Saved values of the mixture self diffusion coefficients
vector_fp m_selfDiffMix;
//! Saved species index of the mixture correlated to self diffusion coefficients
std::vector<std::string> m_selfDiffMixIndex;
//! work space
/*!
* Length is equal to m_nsp
@ -1115,6 +1341,46 @@ namespace Cantera {
//! are current wrt the concentration
bool m_visc_conc_ok;
//! Boolean indicating that the top-level mixture ionic conductivity is current
/*!
* This is turned false for every change in T, P, or C.
*/
bool m_ionCond_mix_ok;
//! Boolean indicating that weight factors wrt ionic conductivty is current
bool m_ionCond_temp_ok;
//! Flag to indicate that the pure species ionic conductivities
//! are current wrt the concentration
bool m_ionCond_conc_ok;
bool m_cond_mix_ok;
//! Boolean indicating that the top-level mixture mobility ratio is current
/*!
* This is turned false for every change in T, P, or C.
*/
bool m_mobRat_mix_ok;
//! Boolean indicating that weight factors wrt mobility ratio is current
bool m_mobRat_temp_ok;
//! Flag to indicate that the pure species mobility ratios
//! are current wrt the concentration
bool m_mobRat_conc_ok;
//! Boolean indicating that the top-level mixture self diffusion is current
/*!
* This is turned false for every change in T, P, or C.
*/
bool m_selfDiff_mix_ok;
//! Boolean indicating that weight factors wrt self diffusion is current
bool m_selfDiff_temp_ok;
//! Flag to indicate that the pure species self diffusion
//! are current wrt the concentration
bool m_selfDiff_conc_ok;
//! Boolean indicating that mixture diffusion coeffs are current
bool m_radi_mix_ok;
@ -1134,10 +1400,10 @@ namespace Cantera {
//! Flag to indicate that the pure species conductivities
//! are current wrt the temperature
bool m_cond_temp_ok;
bool m_lambda_temp_ok;
//! Boolean indicating that mixture conductivity is current
bool m_cond_mix_ok;
bool m_lambda_mix_ok;
//! Mode indicator for transport models -- currently unused.
int m_mode;

100
Cantera/src/transport/LiquidTransportData.cpp
View file

@ -10,6 +10,7 @@
*/
#include "LiquidTransportData.h"
using namespace std;
namespace Cantera {
@ -74,6 +75,11 @@ namespace Cantera {
speciesName = right.speciesName;
hydroRadius = right.hydroRadius;
viscosity = right.viscosity;
ionConductivity = right.ionConductivity;
mobilityRatio = right.mobilityRatio;
mobRatIndex = right.mobRatIndex;
selfDiffusion = right.selfDiffusion;
selfDiffIndex = right.selfDiffIndex;
thermalCond = right.thermalCond;
electCond = right.electCond;
speciesDiffusivity = right.speciesDiffusivity;
@ -167,7 +173,9 @@ namespace Cantera {
thermo_t* thermo ) :
LTPspecies( propNode, name, tp_ind, thermo)
{
m_model = LTR_MODEL_ARRHENIUS;
m_model = LTR_MODEL_ARRHENIUS;
m_temp = 0.0;
m_prop = 0.0;
doublereal A_k, n_k, Tact_k;
getArrhenius(propNode, A_k, n_k, Tact_k);
@ -268,14 +276,16 @@ namespace Cantera {
thermo_t* thermo ) :
LTPspecies( propNode, name, tp_ind, thermo)
{
m_model = LTR_MODEL_POLY;
m_model = LTR_MODEL_POLY;
m_temp = 0.0;
m_prop = 0.0;
getFloatArray(propNode, m_coeffs, true); // if units labeled, convert Angstroms -> meters
getFloatArray(propNode, m_coeffs, "true", "toSI");
if (m_coeffs[0] <= 0.0) {
/* if (m_coeffs[0] <= 0.0) {
throw LTPError("negative or zero " + propNode.name() );
}
}*/
}
//! Copy constructor
@ -308,13 +318,89 @@ namespace Cantera {
doublereal LTPspecies_Poly::getSpeciesTransProp( ) {
doublereal t = m_thermo->temperature();
if (t != m_temp) {
double tempN = 1.0;
if (t != m_temp) {
m_temp=t;
double tempN = 1.0;
for (int i = 0; i < (int) m_coeffs.size() ; i++) {
m_prop += m_coeffs[i] * tempN;
//cout << "m_coeff = " <<m_coeffs[i] << ", tempN = " << tempN << ", m_prop = " << m_prop << endl;
tempN *= m_temp;
}
}
//cout << "m_prop = " << m_prop << endl;
return m_prop;
}
///////////////////////////////////////////////////////////////
//! Construct an LTPspecies object for a liquid tranport property
//! expressed as an exponential in temperature.
/** The transport property is constructed from the XML node,
* \verbatim <propNode>, \endverbatim that is a child of the
* \verbatim <transport> \endverbatim node and specifies a type of
* transport property (like viscosity)
*/
LTPspecies_ExpT::LTPspecies_ExpT( const XML_Node &propNode,
std::string name,
TransportPropertyList tp_ind,
thermo_t* thermo ) :
LTPspecies( propNode, name, tp_ind, thermo)
{
m_model = LTR_MODEL_EXPT;
m_temp = 0.0;
m_prop = 0.0;
getFloatArray(propNode, m_coeffs, "true", "toSI");
/* if (m_coeffs[0] <= 0.0) {
throw LTPError("negative or zero " + propNode.name() );
}*/
}
//! Copy constructor
LTPspecies_ExpT::LTPspecies_ExpT( const LTPspecies_ExpT &right )
: LTPspecies()
{
*this = right; //use assignment operator to do other work
}
//! Assignment operator
LTPspecies_ExpT& LTPspecies_ExpT::operator=(const LTPspecies_ExpT& right )
{
if (&right != this) {
//LTPspecies::operator=(right);
m_speciesName = right.m_speciesName;
m_property = right.m_property;
m_model = right.m_model;
m_coeffs = right.m_coeffs;
m_thermo = right.m_thermo;
m_mixWeight = right.m_mixWeight;
m_temp = right.m_temp;
m_prop = right.m_prop;
}
return *this;
}
//! Return the value for this transport property evaluated
//! from the exponential in temperature expression
doublereal LTPspecies_ExpT::getSpeciesTransProp( ) {
doublereal t = m_thermo->temperature();
if (t != m_temp) {
m_temp=t;
m_prop=m_coeffs[0];
double tempN = 1.0;
for (int i = 1; i < (int) m_coeffs.size() ; i++) {
tempN *= m_temp;
m_prop *= exp(m_coeffs[i] * tempN);
//cout << "m_coeff = " <<m_coeffs[i] << ", tempN = " << tempN << ", m_prop = " << m_prop << endl;
}
}
//cout << "m_prop = " << m_prop << endl;
return m_prop;
}

78
Cantera/src/transport/LiquidTransportData.h
View file

@ -45,6 +45,9 @@ namespace Cantera {
enum TransportPropertyList {
TP_UNKNOWN = -1,
TP_VISCOSITY = 0,
TP_IONCONDUCTIVITY,
TP_MOBILITYRATIO,
TP_SELFDIFFUSION,
TP_THERMALCOND,
TP_DIFFUSIVITY,
TP_HYDRORADIUS,
@ -62,7 +65,8 @@ namespace Cantera {
LTR_MODEL_NOTSET=-1,
LTR_MODEL_CONSTANT,
LTR_MODEL_ARRHENIUS,
LTR_MODEL_POLY
LTR_MODEL_POLY,
LTR_MODEL_EXPT
};
@ -196,6 +200,17 @@ namespace Cantera {
//! Model type for the viscosity
LTPspecies* viscosity;
//! Model type for the ionic conductivity
LTPspecies* ionConductivity;
//! Model type for the mobility ratio
std::vector<LTPspecies*> mobilityRatio;
std::vector<std::string> mobRatIndex;
//! Model type for the self diffusion coefficients
std::vector<LTPspecies*> selfDiffusion;
std::vector<std::string> selfDiffIndex;
//! Model type for the thermal conductivity
LTPspecies* thermalCond;
@ -414,6 +429,67 @@ namespace Cantera {
};
//! Class LTPspecies_ExpT holds transport parameters for a
//! specific liquid-phase species when the transport property
//! is expressed as a exponential in temperature.
/**
* As an example of the input required for LTPspecies_ExpT
* consider the following XML fragment
*
* \verbatim
* <species>
* <!-- thermodynamic properties -->
* <transport>
* <thermalConductivity model="expT">
* <floatArray size="2"> 0.6, -15.0e-5 </floatArray>
* </thermalConductivity>
* <!-- other tranport properties -->
* </transport>
* </species>
* \endverbatim
*/
class LTPspecies_ExpT : public LTPspecies{
public:
LTPspecies_ExpT( const XML_Node &propNode,
std::string name,
TransportPropertyList tp_ind,
thermo_t* thermo );
//! Copy constructor
LTPspecies_ExpT( const LTPspecies_ExpT &right );
//! Assignment operator
LTPspecies_ExpT& operator=(const LTPspecies_ExpT& right );
virtual ~LTPspecies_ExpT( ) { }
//! Returns the pure species tranport property
/*!
* The pure species transport property (i.e. pure species viscosity)
* is returned. Any temperature and composition dependence will be
* adjusted internally according to the information provided.
*/
doublereal getSpeciesTransProp( );
protected:
//! temperature from thermo object
doublereal m_temp;
//! most recent evaluation of transport property
doublereal m_prop;
//! Internal model to adjust species-specific properties for composition.
/** Currently just a place holder, but this method could take
* the composition from the thermo object and adjust coefficients
* accoding to some unspecified model.
*/
void adjustCoeffsForComposition( ){ }
};
}
#endif

689
Cantera/src/transport/LiquidTransportParams.cpp
View file

@ -10,6 +10,11 @@
*/
#include "LiquidTransportParams.h"
#include <iostream>
#include "IonsFromNeutralVPSSTP.h"
#include "MargulesVPSSTP.h"
#include <stdlib.h>
using namespace std;
namespace Cantera {
@ -49,19 +54,44 @@ namespace Cantera {
{
}
LiquidTranInteraction::~LiquidTranInteraction(){
int kmax = m_Aij.size();
for ( int k = 0; k < kmax; k++)
if ( m_Aij[k] ) delete m_Aij[k];
kmax = m_Bij.size();
for ( int k = 0; k < kmax; k++)
if ( m_Bij[k] ) delete m_Bij[k];
kmax = m_Hij.size();
for ( int k = 0; k < kmax; k++)
if ( m_Hij[k] ) delete m_Hij[k];
kmax = m_Sij.size();
for ( int k = 0; k < kmax; k++)
if ( m_Sij[k] ) delete m_Sij[k];
}
void LiquidTranInteraction::init( const XML_Node &compModelNode,
thermo_t* thermo )
{
doublereal poly0;
m_thermo = thermo;
int nsp = thermo->nSpecies();
m_Aij.resize( nsp, nsp, 0.0 );
m_Dij.resize( nsp, nsp, 0.0 );
m_Eij.resize( nsp, nsp, 0.0 );
m_Sij.resize( nsp, nsp, 0.0 );
/*
m_Aij.resize( nsp);
m_Bij.resize( nsp);
m_Hij.resize( nsp);
m_Sij.resize( nsp);
for ( int k = 0; k < nsp; k++ ){
(*m_Aij[k]).resize( nsp, nsp, 0.0);
(*m_Bij[k]).resize( nsp, nsp, 0.0);
(*m_Hij[k]).resize( nsp, nsp, 0.0);
(*m_Sij[k]).resize( nsp, nsp, 0.0);
}
*/
std::string speciesA;
std::string speciesB;
@ -83,10 +113,10 @@ namespace Cantera {
if ( jSpecies < 0 )
throw CanteraError("TransportFactory::getLiquidInteractionsTransportData",
"Unknown species " + speciesB );
if ( xmlChild.hasChild( "Aij" ) ) {
/* if ( xmlChild.hasChild( "Aij" ) ) {
m_Aij(iSpecies,jSpecies) = getFloat( xmlChild, "Aij", "toSI" );
m_Aij(jSpecies,iSpecies) = m_Aij(iSpecies,jSpecies) ;
}
}*/
if ( xmlChild.hasChild( "Eij" ) ) {
m_Eij(iSpecies,jSpecies) = getFloat( xmlChild, "Eij", "actEnergy" );
@ -94,11 +124,73 @@ namespace Cantera {
m_Eij(jSpecies,iSpecies) = m_Eij(iSpecies,jSpecies) ;
}
if ( xmlChild.hasChild( "Sij" ) ) {
if ( xmlChild.hasChild( "Aij" ) ) {
vector_fp poly;
poly0 = getFloat( poly, xmlChild, "Aij", "toSI" );
if ( !poly.size() ) poly.push_back(poly0);
while (m_Aij.size()<poly.size()){
DenseMatrix * aTemp = new DenseMatrix();
aTemp->resize( nsp, nsp, 0.0);
m_Aij.push_back(aTemp);
}
for( int i=0; i<(int)poly.size(); i++ ){
(*m_Aij[i])(iSpecies,jSpecies) = poly[i];
//(*m_Aij[i])(jSpecies,iSpecies) = (*m_Aij[i])(iSpecies,jSpecies) ;
}
}
if ( xmlChild.hasChild( "Bij" ) ) {
vector_fp poly;
poly0 = getFloat( poly, xmlChild, "Bij", "toSI" );
if ( !poly.size() ) poly.push_back(poly0);
while (m_Bij.size()<poly.size()){
DenseMatrix * bTemp = new DenseMatrix();
bTemp->resize( nsp, nsp, 0.0);
m_Bij.push_back(bTemp);
}
for( int i=0; i<(int)poly.size(); i++ ){
(*m_Bij[i])(iSpecies,jSpecies) = poly[i];
//(*m_Bij[i])(jSpecies,iSpecies) = (*m_Bij[i])(iSpecies,jSpecies) ;
}
}
if ( xmlChild.hasChild( "Hij" ) ) {
vector_fp poly;
poly0 = getFloat( poly, xmlChild, "Hij", "actEnergy" );
if ( !poly.size() ) poly.push_back(poly0);
while (m_Hij.size()<poly.size()){
DenseMatrix * hTemp = new DenseMatrix();
hTemp->resize( nsp, nsp, 0.0);
m_Hij.push_back(hTemp);
}
for( int i=0; i<(int)poly.size(); i++ ){
(*m_Hij[i])(iSpecies,jSpecies) = poly[i];
(*m_Hij[i])(iSpecies,jSpecies) /= GasConstant;
//(*m_Hij[i])(jSpecies,iSpecies) = (*m_Hij[i])(iSpecies,jSpecies) ;
}
}
if ( xmlChild.hasChild( "Sij" ) ) {
vector_fp poly;
poly0 = getFloat( poly, xmlChild, "Sij", "actEnergy" );
if ( !poly.size() ) poly.push_back(poly0);
while (m_Sij.size()<poly.size()){
DenseMatrix * sTemp = new DenseMatrix();
sTemp->resize( nsp, nsp, 0.0);
m_Sij.push_back(sTemp);
}
for( int i=0; i<(int)poly.size(); i++ ){
(*m_Sij[i])(iSpecies,jSpecies) = poly[i];
(*m_Sij[i])(iSpecies,jSpecies) /= GasConstant;
//(*m_Sij[i])(jSpecies,iSpecies) = (*m_Sij[i])(iSpecies,jSpecies) ;
}
}
/*0 if ( xmlChild.hasChild( "Sij" ) ) {
m_Sij(iSpecies,jSpecies) = getFloat( xmlChild, "Sij", "toSI" );
m_Sij(iSpecies,jSpecies) /= GasConstant;
m_Sij(jSpecies,iSpecies) = m_Sij(iSpecies,jSpecies) ;
}
//m_Sij(jSpecies,iSpecies) = m_Sij(iSpecies,jSpecies) ;
}*/
if ( xmlChild.hasChild( "Dij" ) ) {
m_Dij(iSpecies,jSpecies) = getFloat( xmlChild, "Dij", "toSI" );
@ -121,7 +213,9 @@ namespace Cantera {
m_thermo = right.m_thermo;
//m_trParam = right.m_trParam;
m_Aij = right.m_Aij;
m_Bij = right.m_Bij;
m_Eij = right.m_Eij;
m_Hij = right.m_Hij;
m_Sij = right.m_Sij;
m_Dij = right.m_Dij;
}
@ -173,6 +267,7 @@ namespace Cantera {
doublereal LTI_Solvent::getMixTransProp( doublereal *speciesValues, doublereal *speciesWeight ) {
int nsp = m_thermo->nSpecies();
doublereal temp = m_thermo->temperature();
doublereal molefracs[nsp];
m_thermo->getMoleFractions(molefracs);
@ -181,23 +276,32 @@ namespace Cantera {
//if weightings are specified, use those
if ( speciesWeight ) {
for ( int k = 0; k < nsp; k++) {
//presume that the weighting is set to 1.0 for solvent and 0.0 for everything else.
value += speciesValues[k] * speciesWeight[k];
molefracs[k] = molefracs[k];
// should be: molefracs[k] = molefracs[k]*speciesWeight[k]; for consistency, but weight(solvent)=1?
}
}
else {
//This does not follow directly a solvent model
throw CanteraError("LTI_Solvent::getMixTransProp","You should be specifying the speciesWeight");
/* //This does not follow directly a solvent model
//although if the solvent mole fraction is dominant
//and the other species values are given or zero,
//it should work.
for ( int k = 0; k < nsp; k++) {
value += speciesValues[k] * molefracs[k];
}
}*/
}
for ( int i = 0; i < nsp; i++ )
for ( int j = 0; j < i; j++ )
value += molefracs[i] * molefracs[j] * m_Aij(i,j) ;
for ( int i = 0; i < nsp; i++ ){
//presume that the weighting is set to 1.0 for solvent and 0.0 for everything else.
value += speciesValues[i] * speciesWeight[i];
for ( int j = 0; j < nsp; j++ ){
for ( int k = 0; k < (int)m_Aij.size(); k++ ){
value += molefracs[i]*molefracs[j]*(*m_Aij[k])(i,j)*pow(molefracs[i],k);
}
for ( int k = 0; k < (int)m_Bij.size(); k++ ){
value += molefracs[i]*molefracs[j]*(*m_Bij[k])(i,j)*temp*pow(molefracs[i],k);
} }
}
return value;
}
@ -206,19 +310,29 @@ namespace Cantera {
doublereal LTI_Solvent::getMixTransProp( std::vector<LTPspecies*> LTPptrs ) {
int nsp = m_thermo->nSpecies();
doublereal temp = m_thermo->temperature();
doublereal molefracs[nsp];
m_thermo->getMoleFractions(molefracs);
doublereal value = 0;
for ( int k = 0; k < nsp; k++) {
//presume that the weighting is set to 1.0 for solvent and 0.0 for everything else.
value += LTPptrs[k]->getSpeciesTransProp() * LTPptrs[k]->getMixWeight( ) ;
molefracs[k] = molefracs[k];
// should be: molefracs[k] = molefracs[k]*LTPptrs[k]->getMixWeight(); for consistency, but weight(solvent)=1?
}
for ( int i = 0; i < nsp; i++ )
for ( int j = 0; j < i; j++ )
value += molefracs[i] * molefracs[j] * m_Aij(i,j) ;
for ( int i = 0; i < nsp; i++ ){
//presume that the weighting is set to 1.0 for solvent and 0.0 for everything else.
value += LTPptrs[i]->getSpeciesTransProp() * LTPptrs[i]->getMixWeight();
for ( int j = 0; j < nsp; j++ ){
for ( int k = 0; k < (int)m_Aij.size(); k++ ){
value += molefracs[i]*molefracs[j]*(*m_Aij[k])(i,j)*pow(molefracs[i],k);
}
for ( int k = 0; k < (int)m_Bij.size(); k++ ){
value += molefracs[i]*molefracs[j]*(*m_Bij[k])(i,j)*temp*pow(molefracs[i],k);
}
}
}
return value;
}
@ -229,105 +343,6 @@ namespace Cantera {
doublereal LTI_MoleFracs::getMixTransProp( doublereal *speciesValues, doublereal *speciesWeight ) {
int nsp = m_thermo->nSpecies();
doublereal molefracs[nsp];
m_thermo->getMoleFractions(molefracs);
doublereal value = 0;
//if weightings are specified, use those
if ( speciesWeight ) {
for ( int k = 0; k < nsp; k++) {
value += speciesValues[k] * speciesWeight[k] * molefracs[k];
}
}
else {
for ( int k = 0; k < nsp; k++) {
value += speciesValues[k] * molefracs[k];
}
}
for ( int i = 0; i < nsp; i++ )
for ( int j = 0; j < i; j++ )
value += molefracs[i] * molefracs[j] * m_Aij(i,j) ;
return value;
}
doublereal LTI_MoleFracs::getMixTransProp( std::vector<LTPspecies*> LTPptrs ) {
int nsp = m_thermo->nSpecies();
doublereal molefracs[nsp];
m_thermo->getMoleFractions(molefracs);
doublereal value = 0;
for ( int k = 0; k < nsp; k++) {
value += LTPptrs[k]->getSpeciesTransProp() * LTPptrs[k]->getMixWeight( ) * molefracs[k];
}
for ( int i = 0; i < nsp; i++ )
for ( int j = 0; j < i; j++ )
value += molefracs[i] * molefracs[j] * m_Aij(i,j) ;
return value;
}
doublereal LTI_MassFracs::getMixTransProp(doublereal *speciesValues, doublereal *speciesWeight ) {
int nsp = m_thermo->nSpecies();
doublereal massfracs[nsp];
m_thermo->getMassFractions(massfracs);
doublereal value = 0;
//if weightings are specified, use those
if ( speciesWeight ) {
for ( int k = 0; k < nsp; k++) {
value += speciesValues[k] * speciesWeight[k] * massfracs[k];
}
}
else {
for ( int k = 0; k < nsp; k++) {
value += speciesValues[k] * massfracs[k];
}
}
for ( int i = 0; i < nsp; i++ )
for ( int j = 0; j < i; j++ )
value += massfracs[i] * massfracs[j] * m_Aij(i,j) ;
return value;
}
doublereal LTI_MassFracs::getMixTransProp( std::vector<LTPspecies*> LTPptrs ) {
int nsp = m_thermo->nSpecies();
doublereal massfracs[nsp];
m_thermo->getMassFractions(massfracs);
doublereal value = 0;
for ( int k = 0; k < nsp; k++) {
value += LTPptrs[k]->getSpeciesTransProp() * LTPptrs[k]->getMixWeight( ) * massfracs[k];
}
for ( int i = 0; i < nsp; i++ )
for ( int j = 0; j < i; j++ )
value += massfracs[i] * massfracs[j] * m_Aij(i,j) ;
return value;
}
doublereal LTI_Log_MoleFracs::getMixTransProp( doublereal *speciesValues, doublereal *speciesWeight ) {
int nsp = m_thermo->nSpecies();
doublereal temp = m_thermo->temperature();
doublereal molefracs[nsp];
@ -338,19 +353,157 @@ namespace Cantera {
//if weightings are specified, use those
if ( speciesWeight ) {
for ( int k = 0; k < nsp; k++) {
value += log( speciesValues[k] ) * speciesWeight[k] * molefracs[k];
molefracs[k] = molefracs[k]*speciesWeight[k];
}
}
else {
for ( int k = 0; k < nsp; k++) {
value += log( speciesValues[k] ) * molefracs[k];
throw CanteraError("LTI_MoleFracs::getMixTransProp","You should be specifying the speciesWeight");
}
for ( int i = 0; i < nsp; i++ ){
value += speciesValues[i] * molefracs[i];
for ( int j = 0; j < nsp; j++ ){
for ( int k = 0; k < (int)m_Aij.size(); k++ ){
value += molefracs[i]*molefracs[j]*(*m_Aij[k])(i,j)*pow(molefracs[i],k);
}
for ( int k = 0; k < (int)m_Bij.size(); k++ ){
value += molefracs[i]*molefracs[j]*(*m_Bij[k])(i,j)*temp*pow(molefracs[i],k);
}
}
}
for ( int i = 0; i < nsp; i++ )
for ( int j = 0; j < i; j++ )
value += molefracs[i] * molefracs[j]
* ( m_Sij(i,j) + m_Eij(i,j) / temp );
return value;
}
doublereal LTI_MoleFracs::getMixTransProp( std::vector<LTPspecies*> LTPptrs ) {
int nsp = m_thermo->nSpecies();
doublereal temp = m_thermo->temperature();
doublereal molefracs[nsp];
m_thermo->getMoleFractions(molefracs);
doublereal value = 0;
for ( int k = 0; k < nsp; k++) {
molefracs[k] = molefracs[k]*LTPptrs[k]->getMixWeight();
}
for ( int i = 0; i < nsp; i++ ){
value += LTPptrs[i]->getSpeciesTransProp() * molefracs[i];
for ( int j = 0; j < nsp; j++ ){
for ( int k = 0; k < (int)m_Aij.size(); k++ ){
value += molefracs[i]*molefracs[j]*(*m_Aij[k])(i,j)*pow(molefracs[i],k);
}
for ( int k = 0; k < (int)m_Bij.size(); k++ ){
value += molefracs[i]*molefracs[j]*(*m_Bij[k])(i,j)*temp*pow(molefracs[i],k);
}
}
}
return value;
}
doublereal LTI_MassFracs::getMixTransProp( doublereal *speciesValues, doublereal *speciesWeight ) {
int nsp = m_thermo->nSpecies();
doublereal temp = m_thermo->temperature();
doublereal massfracs[nsp];
m_thermo->getMassFractions(massfracs);
doublereal value = 0;
//if weightings are specified, use those
if ( speciesWeight ) {
for ( int k = 0; k < nsp; k++) {
massfracs[k] = massfracs[k]*speciesWeight[k];
}
}
else {
throw CanteraError("LTI_MassFracs::getMixTransProp","You should be specifying the speciesWeight");
}
for ( int i = 0; i < nsp; i++ ){
value += speciesValues[i] * massfracs[i];
for ( int j = 0; j < nsp; j++ ){
for ( int k = 0; k < (int)m_Aij.size(); k++ ){
value += massfracs[i]*massfracs[j]*(*m_Aij[k])(i,j)*pow(massfracs[i],k);
}
for ( int k = 0; k < (int)m_Bij.size(); k++ ){
value += massfracs[i]*massfracs[j]*(*m_Bij[k])(i,j)*temp*pow(massfracs[i],k);
}
}
}
return value;
}
doublereal LTI_MassFracs::getMixTransProp( std::vector<LTPspecies*> LTPptrs ) {
int nsp = m_thermo->nSpecies();
doublereal temp = m_thermo->temperature();
doublereal massfracs[nsp];
m_thermo->getMassFractions(massfracs);
doublereal value = 0;
for ( int k = 0; k < nsp; k++) {
massfracs[k] = massfracs[k]*LTPptrs[k]->getMixWeight();
}
for ( int i = 0; i < nsp; i++ ){
value += LTPptrs[i]->getSpeciesTransProp() * massfracs[i];
for ( int j = 0; j < nsp; j++ ){
for ( int k = 0; k < (int)m_Aij.size(); k++ ){
value += massfracs[i]*massfracs[j]*(*m_Aij[k])(i,j)*pow(massfracs[i],k);
}
for ( int k = 0; k < (int)m_Bij.size(); k++ ){
value += massfracs[i]*massfracs[j]*(*m_Bij[k])(i,j)*temp*pow(massfracs[i],k);
}
}
}
return value;
}
doublereal LTI_Log_MoleFracs::getMixTransProp( doublereal *speciesValues, doublereal *speciesWeight ) {
int nsp = m_thermo->nSpecies();
doublereal temp = m_thermo->temperature();
doublereal molefracs[nsp];
m_thermo->getMoleFractions(molefracs);
doublereal value = 0;
//if weightings are specified, use those
if ( speciesWeight ) {
for ( int k = 0; k < nsp; k++) {
molefracs[k] = molefracs[k]*speciesWeight[k];
}
}
else{
throw CanteraError("LTI_Log_MoleFracs::getMixTransProp","You probably should have a speciesWeight when you call getMixTransProp to convert ion mole fractions to molecular mole fractions");
}
for ( int i = 0; i < nsp; i++ ){
value += log( speciesValues[i] ) * molefracs[i];
for ( int j = 0; j < nsp; j++ ){
for ( int k = 0; k < (int)m_Hij.size(); k++ ){
value += molefracs[i]*molefracs[j]*(*m_Hij[k])(i,j)/temp*pow(molefracs[i],k);
//cout << "value = " << value << ", m_Sij = " << (*m_Sij[k])(i,j) << ", m_Hij = " << (*m_Hij[k])(i,j) << endl;
}
for ( int k = 0; k < (int)m_Sij.size(); k++ ){
value -= molefracs[i]*molefracs[j]*(*m_Sij[k])(i,j)*pow(molefracs[i],k);
//cout << "value = " << value << ", m_Sij = " << (*m_Sij[k])(i,j) << ", m_Hij = " << (*m_Hij[k])(i,j) << endl;
}
}
}
value = exp( value );
return value;
@ -363,20 +516,36 @@ namespace Cantera {
doublereal temp = m_thermo->temperature();
doublereal molefracs[nsp];
m_thermo->getMoleFractions(molefracs);
doublereal value = 0;
//if weightings are specified, use those
for ( int k = 0; k < nsp; k++) {
value += log( LTPptrs[k]->getSpeciesTransProp() ) * LTPptrs[k]->getMixWeight() * molefracs[k];
molefracs[k] = molefracs[k]*LTPptrs[k]->getMixWeight( );
}
for ( int i = 0; i < nsp; i++ )
for ( int j = 0; j < i; j++ )
value += molefracs[i] * molefracs[j]
* ( m_Sij(i,j) + m_Eij(i,j) / temp );
for ( int i = 0; i < nsp; i++ ){
value += log( LTPptrs[i]->getSpeciesTransProp() ) * molefracs[i];
for ( int j = 0; j < nsp; j++ ){
for ( int k = 0; k < (int)m_Hij.size(); k++ ){
value += molefracs[i]*molefracs[j]*(*m_Hij[k])(i,j)/temp*pow(molefracs[i],k);
//cout << "1 = " << molefracs[i]+molefracs[j] << endl;
//cout << "value = " << value << ", m_Sij = " << (*m_Sij[k])(i,j) << ", m_Hij = " << (*m_Hij[k])(i,j) << endl;
}
for ( int k = 0; k < (int)m_Sij.size(); k++ ){
value -= molefracs[i]*molefracs[j]*(*m_Sij[k])(i,j)*pow(molefracs[i],k);
//cout << "1 = " << molefracs[i]+molefracs[j] << endl;
//cout << "value = " << value << ", m_Sij = " << (*m_Sij[k])(i,j) << ", m_Hij = " << (*m_Hij[k])(i,j) << endl;
}
}
}
value = exp( value );
// cout << ", viscSpeciesA = " << LTPptrs[0]->getSpeciesTransProp() << endl;
//cout << ", viscSpeciesB = " << LTPptrs[1]->getSpeciesTransProp() << endl;
//cout << "value = " << value << " FINAL" << endl;
return value;
}
@ -422,7 +591,7 @@ namespace Cantera {
return value;
}
void LTI_Pairwise_Interaction::getMatrixTransProp(DenseMatrix &mat, doublereal *speciesValues) {
void LTI_Pairwise_Interaction::getMatrixTransProp( DenseMatrix &mat, doublereal *speciesValues ) {
int nsp = m_thermo->nSpecies();
doublereal temp = m_thermo->temperature();
@ -439,6 +608,205 @@ namespace Cantera {
mat(i,i) = m_diagonals[i]->getSpeciesTransProp() ;
}
void LTI_StefanMaxwell_PPN::setParameters( LiquidTransportParams& trParam ) {
int nsp = m_thermo->nSpecies();
int nBinInt = nsp*(nsp-1)/2;
//vector<std
m_ionCondMix = 0;
m_ionCondMixModel = trParam.ionConductivity;
m_ionCondSpecies.resize(nsp,0);
m_mobRatMix.resize(nsp,nsp,0.0);
m_mobRatMixModel.resize(nBinInt);
m_mobRatSpecies.resize(nBinInt);
m_mobRatIndex.resize(nBinInt);
m_selfDiffMix.resize(nsp,0.0);
m_selfDiffMixModel.resize(nsp);
m_selfDiffSpecies.resize(nsp);
m_selfDiffIndex.resize(nsp);
for ( int k = 0; k < nBinInt; k++ ) {
m_mobRatMixModel[k] = trParam.mobilityRatio[k];
m_mobRatSpecies[k].resize(nsp,0);
m_mobRatIndex[k] = trParam.mobRatIndex[k];
}
for ( int k = 0; k < nsp; k++ ) {
m_selfDiffMixModel[k] = trParam.selfDiffusion[k];
m_selfDiffSpecies[k].resize(nsp,0);
m_selfDiffIndex[k] = trParam.selfDiffIndex[k];
}
for (int k = 0; k < nsp; k++) {
Cantera::LiquidTransportData &ltd = trParam.LTData[k];
m_ionCondSpecies[k] = ltd.ionConductivity;
for ( int j = 0; j < nBinInt; j++ ){
m_mobRatSpecies[j][k] = ltd.mobilityRatio[j];
}
for ( int j = 0; j < nsp; j++ ){
m_selfDiffSpecies[j][k] = ltd.selfDiffusion[j];
}
}
}
doublereal LTI_StefanMaxwell_PPN::getMixTransProp( doublereal *speciesValues, doublereal *speciesWeight ) {
int nsp = m_thermo->nSpecies();
doublereal molefracs[nsp];
m_thermo->getMoleFractions(molefracs);
doublereal value = 0;
throw LTPmodelError( "Calling LTI_StefanMaxwell_PPN::getMixTransProp does not make sense." );
return value;
}
doublereal LTI_StefanMaxwell_PPN::getMixTransProp( std::vector<LTPspecies*> LTPptrs ) {
int nsp = m_thermo->nSpecies();
doublereal molefracs[nsp];
m_thermo->getMoleFractions(molefracs);
doublereal value = 0;
throw LTPmodelError( "Calling LTI_StefanMaxwell_PPN::getMixTransProp does not make sense." );
return value;
}
void LTI_StefanMaxwell_PPN::getMatrixTransProp( DenseMatrix &mat, doublereal *speciesValues ) {
//CAL
IonsFromNeutralVPSSTP * ions_thermo = dynamic_cast<IonsFromNeutralVPSSTP *>(m_thermo);
int i, j, k;
int nsp = m_thermo->nSpecies();
if (nsp != 3) throw CanteraError("LTI_StefanMaxwell_PPN::getMatrixTransProp","Function may only be called with a 3-ion system");
int nBinInt = nsp*(nsp-1)/2;
doublereal temp = m_thermo->temperature();
doublereal molefracs[nsp];
m_thermo->getMoleFractions( molefracs );
vector_fp neut_molefracs;
ions_thermo->getNeutralMolecMoleFractions(neut_molefracs);
vector<int> cation;
vector<int> anion;
ions_thermo->getCationList(cation);
ions_thermo->getAnionList(anion);
vector<std::string> speciesNames;
ions_thermo->getSpeciesNames(speciesNames);
// Reaction Coeffs and Charges
std::vector<double> viS(6);
std::vector<double> charges(3);
ions_thermo->getDissociationCoeffs(viS,charges);
if ((int)anion.size() != 1)
throw CanteraError("LTI_StefanMaxwell_PPN::getMatrixTransProp","Must have one anion only for StefanMaxwell_PPN");
if ((int)cation.size() != 2)
throw CanteraError("LTI_StefanMaxwell_PPN::getMatrixTransProp","Must have two cations of equal charge for StefanMaxwell_PPN");
if (charges[cation[0]] != charges[cation[1]])
throw CanteraError("LTI_StefanMaxwell_PPN::getMatrixTransProp","Cations must be of equal charge for StefanMaxwell_PPN");
/*
cout << "cation 0: " << speciesNames[cation[0]] << endl;
cout << "cation 1: " << speciesNames[cation[1]] << endl;
cout << "anion 0: " << speciesNames[anion[0]] << endl;
*/
doublereal kappa = m_ionCondMixModel->getMixTransProp(m_ionCondSpecies);
MargulesVPSSTP * marg_thermo = dynamic_cast<MargulesVPSSTP *> (ions_thermo->neutralMoleculePhase_);
doublereal vol = m_thermo->molarVolume();
typedef std::vector<int> intVec;
std::vector<intVec> mobRatIndexMap;
mobRatIndexMap.resize(nsp);
for ( k = 0; k < nsp; k++ )
mobRatIndexMap[k].resize(nsp,0);
for ( k = 0; k < nBinInt; k++ ) {
bool missingIndex = 1;
for ( i = 0; i < nsp; i++ ) {
for ( j = 0; j < i; j++ ) {
if ( (speciesNames[i] + ":" + speciesNames[j]) == m_mobRatIndex[k] ) {
mobRatIndexMap[i][j] = k;
m_mobRatMix(i,j) = m_mobRatMixModel[k]->getMixTransProp( m_mobRatSpecies[k] );
if ( m_mobRatMix(i,j) > 0 ) m_mobRatMix(j,i) = 1.0/m_mobRatMix(i,j);
missingIndex = 0;
break;
}
else if ( (speciesNames[j] + ":" + speciesNames[i]) == m_mobRatIndex[k] ) {
m_mobRatMix(j,i) = m_mobRatMixModel[k]->getMixTransProp( m_mobRatSpecies[k] );
if ( m_mobRatMix(j,i) > 0 ) m_mobRatMix(i,j) = 1.0/m_mobRatMix(j,i);
missingIndex = 0;
break;
}
}
}
if ( missingIndex ) throw CanteraError("LTI_StefanMaxwell_PPN::getMixTransProp","Incorrect names for mobility ratio of " + m_mobRatIndex[k] + " rather than i.e. " + speciesNames[0] + ":" + speciesNames[1]);
}
for ( k = 0; k < nsp; k++ ){
j = 0;
while (m_selfDiffIndex[k] != speciesNames[j]) {
j++;
if (j == nsp) throw CanteraError("LTI_StefanMaxwell_PPN::getMixTransProp","Incorrect names for self diffusion of " + m_selfDiffIndex[k] + " rather than i.e. " + speciesNames[0]);
}
m_selfDiffMix[j] = m_selfDiffMixModel[k]->getMixTransProp( m_selfDiffSpecies[k] );
}
/*
for ( i = 0; i < nsp; i++ ) {
cout << "D" << i << "* = " << m_selfDiffMix[i] << endl;
for ( j = 0; j < nsp; j++ ) {
cout << "ratio" << i << j << " = " << m_mobRatMix(i,j) << endl;
}
}
*/
int vP = max(viS[cation[0]],viS[cation[1]]);
int vM = viS[anion[0]];
int zP = charges[cation[0]];
int zM = charges[anion[0]];
doublereal xA, xB, eps;
doublereal inv_vP_vM_MutualDiff;
vector_fp dlnActCoeffdlnN;
dlnActCoeffdlnN.resize(neut_molefracs.size(),0.0);
std::string cationIndex (4,'0');
for ( i = 0; i < 2; i++ )
for ( j = 0; j < 2; j++ )
if ( viS[i*nsp+cation[j]] > 0 )
cationIndex[i*2+j] = '1';
if ( (cationIndex == "1001") | (cationIndex == "0110") ) {
xA = neut_molefracs[cation[0]];
xB = neut_molefracs[cation[1]];
eps = (1-m_mobRatMix(cation[1],cation[0]))/(xA+xB*m_mobRatMix(cation[1],cation[0]));
marg_thermo->getdlnActCoeffdlnN(&dlnActCoeffdlnN[0]);
inv_vP_vM_MutualDiff = (xA*(1+dlnActCoeffdlnN[cation[1]])/m_selfDiffMix[cation[1]]+xB*(1+dlnActCoeffdlnN[cation[0]])/m_selfDiffMix[cation[0]]);
//marg_thermo->getdlnActCoeffdlnX(&dlnActCoeffdlnN[0]);
//inv_vP_vM_MutualDiff = (xA*(1+dlnActCoeffdlnN[cation[1]])/m_selfDiffMix[cation[1]]+xB*(1+dlnActCoeffdlnN[cation[0]])/m_selfDiffMix[cation[0]]);
}
else
throw CanteraError("LTI_StefanMaxwell_PPN::getMixTransProp","Dissociation reactions don't make sense: cationIndex = " + cationIndex);
mat.resize( nsp, nsp, 0.0 );
mat(cation[0],cation[1]) = mat(cation[1],cation[0]) = (1+vM/vP)*(1+eps*xB)*(1-eps*xA)*inv_vP_vM_MutualDiff-zP*zP*Faraday*Faraday/GasConstant/temp/kappa/vol;
mat(cation[0],anion[0]) = mat(anion[0],cation[0]) = (1+vP/vM)*(-eps*xB*(1-eps*xA)*inv_vP_vM_MutualDiff)-zP*zM*Faraday*Faraday/GasConstant/temp/kappa/vol;
mat(cation[1],anion[0]) = mat(anion[0],cation[1]) = (1+vP/vM)*(eps*xA*(1+eps*xB)*inv_vP_vM_MutualDiff)-zP*zM*Faraday*Faraday/GasConstant/temp/kappa/vol;
for ( i = 0; i < nsp; i++ ) {
for ( j = 0; j < nsp; j++ ) {
mat(i,j) = 1.0/mat(i,j);
//cout << "D" << i << j << " = " << mat(i,j) << endl;
}
}
}
doublereal LTI_StokesEinstein::getMixTransProp( doublereal *speciesValues, doublereal *speciesWeight ) {
@ -503,6 +871,61 @@ namespace Cantera {
delete viscSpec;
}
doublereal LTI_MoleFracs_ExpT::getMixTransProp( doublereal *speciesValues, doublereal *speciesWeight ) {
int nsp = m_thermo->nSpecies();
doublereal temp = m_thermo->temperature();
doublereal molefracs[nsp];
m_thermo->getMoleFractions(molefracs);
doublereal value = 0;
//if weightings are specified, use those
if ( speciesWeight ) {
for ( int k = 0; k < nsp; k++) {
molefracs[k] = molefracs[k]*speciesWeight[k];
}
}
else {
throw CanteraError("LTI_MoleFracs_ExpT::getMixTransProp","You should be specifying the speciesWeight");
}
for ( int i = 0; i < nsp; i++ ){
value += speciesValues[i] * molefracs[i];
for ( int j = 0; j < nsp; j++ ){
for ( int k = 0; k < (int)m_Aij.size(); k++ ){
value += molefracs[i]*molefracs[j]*(*m_Aij[k])(i,j)*pow(molefracs[i],k)*exp((*m_Bij[k])(i,j)*temp);
}
}
}
return value;
}
doublereal LTI_MoleFracs_ExpT::getMixTransProp( std::vector<LTPspecies*> LTPptrs ) {
int nsp = m_thermo->nSpecies();
doublereal temp = m_thermo->temperature();
doublereal molefracs[nsp];
m_thermo->getMoleFractions(molefracs);
doublereal value = 0;
for ( int k = 0; k < nsp; k++) {
molefracs[k] = molefracs[k]*LTPptrs[k]->getMixWeight();
}
for ( int i = 0; i < nsp; i++ ){
value += LTPptrs[i]->getSpeciesTransProp() * molefracs[i];
for ( int j = 0; j < nsp; j++ ){
for ( int k = 0; k < (int)m_Aij.size(); k++ ){
value += molefracs[i]*molefracs[j]*(*m_Aij[k])(i,j)*pow(molefracs[i],k)*exp((*m_Bij[k])(i,j)*temp);
}
}
}
return value;
}

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

@ -33,57 +33,62 @@ namespace Cantera {
};
//! Composition dependence type for liquid mixture transport properties
/*!
* Types of temperature dependencies:
* - 0 - Mixture calculations with this property are not allowed
* - 1 - Use solvent (species 0) properties
* - 2 - Properties weighted linearly by mole fractions
* - 3 - Properties weighted linearly by mass fractions
* - 4 - Properties weighted logarithmically by mole fractions (interaction energy weighting)
* - 5 - Interactions given pairwise between each possible species (i.e. D_ij)
*
* \verbatim
* <transport model="Liquid">
* <viscosity>
* <compositionDependence model="logMoleFractions">
* <interaction>
* <speciesA> LiCl(L) </speciesA>
* <speciesB> KCl(L) </speciesB>
* <Eij units="J/kmol"> -1.0 </Eij>
* <Sij units="J/kmol/K"> 1.0E-1 </Eij>
* </interaction>
* </compositionDependence>
* </viscosity>
* <speciesDiffusivity>
* <compositionDependence model="pairwiseInteraction">
* <interaction>
* <speciesA> Li+ </speciesA>
* <speciesB> K+ </speciesB>
* <Dij units="m2/s"> 1.5 </Dij>
* </interaction>
* <interaction>
* <speciesA> K+ </speciesA>
* <speciesB> Cl- </speciesB>
* <Dij units="m2/s"> 1.0 </Dij>
* </interaction>
* <interaction>
* <speciesA> Li+ </speciesA>
* <speciesB> Cl- </speciesB>
* <Dij units="m2/s"> 1.2 </Dij>
* </interaction>
* </compositionDependence>
* </speciesDiffusivity>
* <thermalConductivity>
* <compositionDependence model="massFractions"/>
* </thermalConductivity>
* <hydrodynamicRadius>
* <compositionDependence model="none"/>
* </hydrodynamicRadius>
* </transport>
* \endverbatim
*
*/
//! Composition dependence type for liquid mixture transport properties
/*!
* Types of temperature dependencies:
* - 0 - Mixture calculations with this property are not allowed
* - 1 - Use solvent (species 0) properties
* - 2 - Properties weighted linearly by mole fractions
* - 3 - Properties weighted linearly by mass fractions
* - 4 - Properties weighted logarithmically by mole fractions (interaction energy weighting)
* - 5 - Interactions given pairwise between each possible species (i.e. D_ij)
*
* \verbatim
* <transport model="Liquid">
* <viscosity>
* <compositionDependence model="logMoleFractions">
* <interaction>
* <speciesA> LiCl(L) </speciesA>
* <speciesB> KCl(L) </speciesB>
* <Eij units="J/kmol"> -1.0 </Eij>
* <Sij units="J/kmol/K"> 1.0E-1 </Sij>
* -or- <Sij>
* <floatArray units="J/kmol/K"> 1.0E-1, 0.001 0.01 </floatArray>
* </Sij>
* -same form for Hij,Aij,Bij-
* </interaction>
* </compositionDependence>
* </viscosity>
* <speciesDiffusivity>
* <compositionDependence model="pairwiseInteraction">
* <interaction>
* <speciesA> Li+ </speciesA>
* <speciesB> K+ </speciesB>
* <Dij units="m2/s"> 1.5 </Dij>
* </interaction>
* <interaction>
* <speciesA> K+ </speciesA>
* <speciesB> Cl- </speciesB>
* <Dij units="m2/s"> 1.0 </Dij>
* </interaction>
* <interaction>
* <speciesA> Li+ </speciesA>
* <speciesB> Cl- </speciesB>
* <Dij units="m2/s"> 1.2 </Dij>
* </interaction>
* </compositionDependence>
* </speciesDiffusivity>
* <thermalConductivity>
* <compositionDependence model="massFractions"/>
* </thermalConductivity>
* <hydrodynamicRadius>
* <compositionDependence model="none"/>
* </hydrodynamicRadius>
* </transport>
* \endverbatim
*
*/
enum LiquidTranMixingModel {
LTI_MODEL_NOTSET=-1,
LTI_MODEL_NONE,
@ -92,7 +97,9 @@ namespace Cantera {
LTI_MODEL_MASSFRACS,
LTI_MODEL_LOG_MOLEFRACS,
LTI_MODEL_PAIRWISE_INTERACTION,
LTI_MODEL_STOKES_EINSTEIN
LTI_MODEL_STEFANMAXWELL_PPN,
LTI_MODEL_STOKES_EINSTEIN,
LTI_MODEL_MOLEFRACS_EXPT
};
@ -132,7 +139,7 @@ namespace Cantera {
LiquidTranInteraction& operator=( const LiquidTranInteraction &right );
//! destructor
virtual ~LiquidTranInteraction() { ; }
virtual ~LiquidTranInteraction();
//! initialize LiquidTranInteraction objects with thermo and XML node
/**
@ -172,14 +179,21 @@ namespace Cantera {
//LiquidTransportParams* m_trParam;
//! Matrix of interactions (no temperature dependence, dimensionless)
DenseMatrix m_Aij;
//! Matrix of interaction coefficients for polynomial in molefraction*weight of speciesA (no temperature dependence, dimensionless)
std::vector<DenseMatrix*> m_Aij;
//! Matrix of interaction coefficients for polynomial in molefraction*weight of speciesA (linear temperature dependence, units 1/K)
std::vector<DenseMatrix*> m_Bij;
//! Matrix of interactions (in energy units, 1/RT temperature dependence)
DenseMatrix m_Eij;
//! Matrix of interactions (in entropy units, divided by R)
DenseMatrix m_Sij;
//! Matrix of interaction coefficients for polynomial in molefraction*weight of speciesA (in energy units, 1/RT temperature dependence)
std::vector<DenseMatrix*> m_Hij;
//! Matrix of interaction coefficients for polynomial in molefraction*weight of speciesA (in entropy units, divided by R)
std::vector<DenseMatrix*> m_Sij;
//DenseMatrix m_Sij;
//! Matrix of interactions
DenseMatrix m_Dij;
@ -205,6 +219,11 @@ namespace Cantera {
std::vector<Cantera::LiquidTransportData> LTData;
LiquidTranInteraction* viscosity;
LiquidTranInteraction* ionConductivity;
std::vector<LiquidTranInteraction*> mobilityRatio;
std::vector<std::string> mobRatIndex;
std::vector<LiquidTranInteraction*> selfDiffusion;
std::vector<std::string> selfDiffIndex;
LiquidTranInteraction* thermalCond;
LiquidTranInteraction* speciesDiffusivity;
LiquidTranInteraction* electCond;
@ -213,7 +232,19 @@ namespace Cantera {
//! Model for species interaction effects for viscosity
//! Takes enum LiquidTranMixingModel
LiquidTranMixingModel model_viscosity;
//! Model for species interaction effects for ionic conductivity
//! Takes enum LiquidTranMixingModel
LiquidTranMixingModel model_ionConductivity;
//! Model for species interaction effects for mobility ratio
//! Takes enum LiquidTranMixingModel
std::vector<LiquidTranMixingModel*> model_mobilityRatio;
//! Model for species interaction effects for mobility ratio
//! Takes enum LiquidTranMixingModel
std::vector<LiquidTranMixingModel*> model_selfDiffusion;
//! Interaction associated with linear weighting of
//! thermal conductivity.
/**
@ -233,7 +264,7 @@ namespace Cantera {
//! Interaction associated with linear weighting of
//! thermal conductivity.
/**
* This is used for either LTI_MODEL_PAIRWISE_INTERACTION.
* This is used for either LTI_MODEL_PAIRWISE_INTERACTION or LTI_MODEL_STEFANMAXWELL_PPN.
* These provide species interaction coefficients associated with
* the Stefan-Maxwell formulation.
*/
@ -256,7 +287,8 @@ namespace Cantera {
class LTI_MassFracs;
class LTI_Log_MoleFracs;
class LTI_Pairwise_Interaction;
class LTI_StefanMaxwell_PPN;
class LTI_MoleFracs_ExpT;
class LTI_Solvent : public LiquidTranInteraction {
public:
@ -283,7 +315,7 @@ namespace Cantera {
doublereal getMixTransProp( doublereal *valueSpecies, doublereal *weightSpecies = 0 );
doublereal getMixTransProp( std::vector<LTPspecies*> LTPptrs ) ;
void getMatrixTransProp( DenseMatrix &mat, doublereal* speciesValues = 0 ) { mat = m_Aij; }
void getMatrixTransProp( DenseMatrix &mat, doublereal* speciesValues = 0 ) { mat = (*m_Aij[0]); }
protected:
@ -325,7 +357,7 @@ namespace Cantera {
doublereal getMixTransProp( doublereal *valueSpecies, doublereal *weightSpecies = 0 );
doublereal getMixTransProp( std::vector<LTPspecies*> LTPptrs ) ;
void getMatrixTransProp( DenseMatrix &mat, doublereal* speciesValues = 0 ) { mat = m_Aij; }
void getMatrixTransProp( DenseMatrix &mat, doublereal* speciesValues = 0 ) { mat = (*m_Aij[0]); }
protected:
@ -368,7 +400,7 @@ namespace Cantera {
doublereal getMixTransProp( doublereal *valueSpecies, doublereal *weightSpecies = 0 );
doublereal getMixTransProp( std::vector<LTPspecies*> LTPptrs ) ;
void getMatrixTransProp( DenseMatrix &mat, doublereal* speciesValues = 0 ) { mat = m_Aij; }
void getMatrixTransProp( DenseMatrix &mat, doublereal* speciesValues = 0 ) { mat = (*m_Aij[0]); }
protected:
@ -502,13 +534,84 @@ namespace Cantera {
doublereal getMixTransProp( std::vector<LTPspecies*> LTPptrs ) ;
void getMatrixTransProp( DenseMatrix &mat, doublereal* speciesValues = 0 ) ;
protected:
std::vector<LTPspecies*> m_diagonals;
};
//! Stefan Maxwell Diffusion Coefficients can be solved for given
//! ion conductivity, mobility ratios, and self diffusion coeffs.
//! This method is only valid for a common anion mixture of two
//! salts with cations of equal charge.
/**
* This class requres you specify
* 1 - ion conductivity
* 2 - mobility ratio of the two cations
* 3 - mutual diffusion coefficient (can be approximated using
* the self diffusion coefficients of the cations
*
* Sample input for this method is
* \verbatim
* <transport model="Liquid">
* <speciesDiffusivity>
* <compositionDependence model="stefanMaxwell_PPN">
* </compositionDependence>
* </speciesDiffusivity>
* </transport>
* \endverbatim
*
*/
class LTI_StefanMaxwell_PPN : public LiquidTranInteraction {
public:
LTI_StefanMaxwell_PPN( TransportPropertyList tp_ind = TP_UNKNOWN ) :
LiquidTranInteraction( tp_ind )
{
m_model = LTI_MODEL_STEFANMAXWELL_PPN;
}
//! Copy constructor
// LTI_StefanMaxwell_PPN( const LTI_StefanMaxwell_PPN &right );
//! Assignment operator
// LTI_StefanMaxwell_PPN& operator=( const LTI_StefanMaxwell_PPN &right );
virtual ~LTI_StefanMaxwell_PPN( ) { }
void setParameters( LiquidTransportParams& trParam ) ;
//! Return the mixture transport property value.
/**
* Takes the separate species transport properties
* as input (this method does not know what
* transport property it is at this point.
*/
doublereal getMixTransProp( doublereal *valueSpecies, doublereal *weightSpecies = 0 );
doublereal getMixTransProp( std::vector<LTPspecies*> LTPptrs ) ;
void getMatrixTransProp( DenseMatrix &mat, doublereal* speciesValues = 0 ) ;
//CAL void getMatrixTransProp( DenseMatrix &mat, LiquidTransport* lt, doublereal* speciesValues = 0 ) ;
protected:
doublereal m_ionCondMix;
LiquidTranInteraction * m_ionCondMixModel;
std::vector<LTPspecies*> m_ionCondSpecies;
typedef std::vector<LTPspecies*> LTPvector;
DenseMatrix m_mobRatMix;
std::vector<LiquidTranInteraction*> m_mobRatMixModel;
std::vector<LTPvector> m_mobRatSpecies;
std::vector<std::string> m_mobRatIndex;
std::vector<LiquidTranInteraction*> m_selfDiffMixModel;
vector_fp m_selfDiffMix;
std::vector<LTPvector> m_selfDiffSpecies;
std::vector<std::string> m_selfDiffIndex;
};
class LTI_StokesEinstein : public LiquidTranInteraction {
public:
@ -539,7 +642,6 @@ namespace Cantera {
doublereal getMixTransProp( std::vector<LTPspecies*> LTPptrs ) ;
void getMatrixTransProp( DenseMatrix &mat, doublereal* speciesValues = 0 ) ;
protected:
std::vector<LTPspecies*> m_viscosity;
@ -547,6 +649,48 @@ namespace Cantera {
};
//! Simple mole fraction weighting of transport properties
/**
* This model weights the transport property by the mole
* fractions.
* The overall formula for the mixture viscosity is
*
* \f[ \eta_{mix} = \sum_i X_i \eta_i
* + \sum_i \sum_j X_i X_j A_{i,j} \f].
*/
class LTI_MoleFracs_ExpT : public LiquidTranInteraction {
public:
LTI_MoleFracs_ExpT( TransportPropertyList tp_ind = TP_UNKNOWN ) :
LiquidTranInteraction( tp_ind )
{
m_model = LTI_MODEL_MOLEFRACS_EXPT;
}
//! Copy constructor
// LTI_MoleFracs_ExpT( const LTI_MoleFracs_ExpT &right );
//! Assignment operator
// LTI_MoleFracs_ExpT& operator=( const LTI_MoleFracs_ExpT &right );
virtual ~LTI_MoleFracs_ExpT( ) { }
//! Return the mixture transport property value.
/**
* Takes the separate species transport properties
* as input (this method does not know what
* transport property it is at this point.
*/
doublereal getMixTransProp( doublereal *valueSpecies, doublereal *weightSpecies = 0 );
doublereal getMixTransProp( std::vector<LTPspecies*> LTPptrs ) ;
void getMatrixTransProp( DenseMatrix &mat, doublereal* speciesValues = 0 ) { mat = (*m_Aij[0]); }
//CAL void getMatrixTransProp( DenseMatrix &mat, LiquidTransport* lt, doublereal* speciesValues = 0 ) { mat = (*m_Aij[0]); }
protected:
};
}

70
Cantera/src/transport/TransportBase.h
View file

@ -28,7 +28,7 @@
#define CT_TRANSPORTBASE_H
#include "ThermoPhase.h"
#include "DenseMatrix.h"
namespace Cantera {
@ -80,6 +80,8 @@ namespace Cantera {
* - VB_MOLEAVG Diffusion velocities are based on the mole averaged velocities
* - VB_SPECIES_0 Diffusion velocities are based on the relative motion wrt species 0
* - VB_SPECIES_1 Diffusion velocities are based on the relative motion wrt species 1
* - VB_SPECIES_2 Diffusion velocities are based on the relative motion wrt species 2
* - VB_SPECIES_3 Diffusion velocities are based on the relative motion wrt species 3
*
* @ingroup tranprops
*/
@ -97,6 +99,10 @@ namespace Cantera {
const VelocityBasis VB_SPECIES_0 = 0;
//! Diffusion velocities are based on the relative motion wrt species 1
const VelocityBasis VB_SPECIES_1 = 1;
//! Diffusion velocities are based on the relative motion wrt species 2
const VelocityBasis VB_SPECIES_2 = 2;
//! Diffusion velocities are based on the relative motion wrt species 3
const VelocityBasis VB_SPECIES_3 = 3;
//@}
//! Base class for transport property managers.
@ -308,6 +314,68 @@ namespace Cantera {
virtual doublereal bulkViscosity()
{ return err("bulkViscosity"); }
/**
* The ionic conducitivity in 1/ohm/m.
*/
virtual doublereal ionConductivity()
{ return err("ionConductivity"); }
//! Returns the pure species ionic conducitivity
/*!
* The units are 1/ohm/m and the length is the number of species
*
* @param ionCond Vector of ionic conductivities
*/
virtual void getSpecoesIonConductivity(doublereal* const ionCond)
{ err("getSpeciesIonConductivity"); }
/**
* The mobility ratio.
*/
virtual void mobilityRatio(vector_fp& mobRat, std::vector<std::string>& mobRatIndex)
{ err("mobilityRatio"); }
virtual void mobilityRatio(double* mobRat, std::vector<std::string>& mobRatIndex)
{ err("mobilityRatio"); }
//! Returns the pure species limit of the mobility ratios
/*!
* The value is dimensionless and the length is the number of species
*
* @param mobRat Vector of mobility ratios
*/
virtual void getSpeciesMobilityRatio(DenseMatrix& mobRat, std::vector<std::string>& mobRatIndex)
{ err("getSpeciesMobilityRatio"); }
virtual void getSpeciesMobilityRatio(double** mobRat, std::vector<std::string>& mobRatIndex)
{ err("getSpeciesMobilityRatio"); }
//! Returns the self diffusion coefficients in the solution
/*!
* The self diffusion calculation is handled by subclasses of
* LiquidTranInteraction as specified in the input file.
* These in turn employ subclasses of LTPspecies to
* determine the individual species self diffusion coeffs.
*/
virtual void selfDiffusion(vector_fp& selfDiff, std::vector<std::string>& selfDiffIndex)
{ err("selfDiffusion"); }
virtual void selfDiffusion(double* selfDiff, std::vector<std::string>& selfDiffIndex)
{ err("selfdiffusion"); }
//! Returns the pure species self diffusion in solution of each species
/*!
* The pure species molar volumes are evaluated using the
* appropriate subclasses of LTPspecies as specified in the
* input file.
*
* @param selfDiff array of length "number of species"
* to hold returned self diffusion coeffs.
*/
virtual void getSpeciesSelfDiffusion(DenseMatrix& selfDiff, std::vector<std::string>& selfDiffIndex)
{ err("getSpeciesSelfDiffusion"); }
virtual void getSpeciesSelfDiffusion(double** selfDiff, std::vector<std::string>& selfDiffIndex)
{ err("getSpeciesSelfDiffusion"); }
//! Returns the mixture thermal conductivity in W/m/K.
/*!
* Units are in W / m K or equivalently kg m / s3 K

138
Cantera/src/transport/TransportFactory.cpp
View file

@ -232,6 +232,9 @@ namespace Cantera {
//m_models["Radiative"] = cRadiative;
m_tranPropMap["viscostiy"] = TP_VISCOSITY;
m_tranPropMap["ionConductivity"] = TP_IONCONDUCTIVITY;
m_tranPropMap["mobilityRatio"] = TP_MOBILITYRATIO;
m_tranPropMap["selfDiffusion"] = TP_SELFDIFFUSION;
m_tranPropMap["thermalConductivity"] = TP_THERMALCOND;
m_tranPropMap["speciesDiffusivity"] = TP_DIFFUSIVITY;
m_tranPropMap["hydrodynamicRadius"] = TP_HYDRORADIUS;
@ -241,6 +244,7 @@ namespace Cantera {
m_LTRmodelMap["constant"] = LTR_MODEL_CONSTANT;
m_LTRmodelMap["arrhenius"] = LTR_MODEL_ARRHENIUS;
m_LTRmodelMap["coeffs"] = LTR_MODEL_POLY;
m_LTRmodelMap["exptemp"] = LTR_MODEL_EXPT;
m_LTImodelMap[""] = LTI_MODEL_NOTSET;
m_LTImodelMap["none"] = LTI_MODEL_NONE;
@ -249,6 +253,8 @@ namespace Cantera {
m_LTImodelMap["massFractions"] = LTI_MODEL_MASSFRACS;
m_LTImodelMap["logMoleFractions"] = LTI_MODEL_LOG_MOLEFRACS;
m_LTImodelMap["pairwiseInteraction"] = LTI_MODEL_PAIRWISE_INTERACTION;
m_LTImodelMap["stefanMaxwell_PPN"] = LTI_MODEL_STEFANMAXWELL_PPN;
m_LTImodelMap["moleFractionsExpT"] = LTI_MODEL_MOLEFRACS_EXPT;
}
/**
@ -312,6 +318,12 @@ namespace Cantera {
tp_ind,
thermo );
break;
case LTR_MODEL_EXPT:
ltps = new LTPspecies_ExpT( trNode,
name,
tp_ind,
thermo );
break;
default:
throw CanteraError("newLTP","unknown transport model: " + model );
ltps = new LTPspecies( trNode,
@ -359,11 +371,20 @@ namespace Cantera {
lti->init( trNode, thermo );
lti->setParameters( trParam );
break;
case LTI_MODEL_STEFANMAXWELL_PPN:
lti = new LTI_StefanMaxwell_PPN( tp_ind );
lti->init( trNode, thermo );
lti->setParameters( trParam );
break;
case LTI_MODEL_STOKES_EINSTEIN:
lti = new LTI_StokesEinstein( tp_ind );
lti->init( trNode, thermo );
lti->setParameters( trParam );
break;
case LTI_MODEL_MOLEFRACS_EXPT:
lti = new LTI_MoleFracs_ExpT( tp_ind );
lti->init( trNode, thermo );
break;
default:
// throw CanteraError("newLTI","unknown transport model: " + model );
lti = new LiquidTranInteraction( tp_ind );
@ -966,7 +987,7 @@ namespace Cantera {
for (int iChild = 0; iChild < num; iChild++) {
XML_Node &xmlChild = trNode.child(iChild);
std::string nodeName = xmlChild.name();
switch ( m_tranPropMap[nodeName] ) {
case TP_VISCOSITY:
data.viscosity = newLTP( xmlChild,
@ -974,23 +995,63 @@ namespace Cantera {
m_tranPropMap[nodeName],
trParam.thermo );
break;
case TP_IONCONDUCTIVITY:
data.ionConductivity = newLTP( xmlChild,
name,
m_tranPropMap[nodeName],
trParam.thermo );
break;
case TP_MOBILITYRATIO:
{
int iSpec;
int numSpec = xmlChild.nChildren();
data.mobRatIndex.resize(numSpec);
data.mobilityRatio.resize(numSpec);
for (iSpec = 0; iSpec< numSpec; iSpec++){
XML_Node &propSpecNode = xmlChild.child(iSpec);
std::string specName = propSpecNode.name();
data.mobRatIndex[iSpec] = specName;
data.mobilityRatio[iSpec] = newLTP( propSpecNode,
name,
m_tranPropMap[nodeName],
trParam.thermo );
};
};
break;
case TP_SELFDIFFUSION:
{
int iSpec;
int numSpec = xmlChild.nChildren();
data.selfDiffIndex.resize(numSpec);
data.selfDiffusion.resize(numSpec);
for (iSpec = 0; iSpec< numSpec; iSpec++){
XML_Node &propSpecNode = xmlChild.child(iSpec);
std::string specName = propSpecNode.name();
data.selfDiffIndex[iSpec] = specName;
data.selfDiffusion[iSpec] = newLTP( propSpecNode,
name,
m_tranPropMap[nodeName],
trParam.thermo );
};
};
break;
case TP_THERMALCOND:
data.thermalCond = newLTP( xmlChild,
name,
m_tranPropMap[nodeName],
trParam.thermo );
name,
m_tranPropMap[nodeName],
trParam.thermo );
break;
case TP_DIFFUSIVITY:
data.speciesDiffusivity = newLTP( xmlChild,
name,
m_tranPropMap[nodeName],
trParam.thermo );
name,
m_tranPropMap[nodeName],
trParam.thermo );
break;
case TP_HYDRORADIUS:
data.hydroRadius = newLTP( xmlChild,
name,
m_tranPropMap[nodeName],
trParam.thermo );
name,
m_tranPropMap[nodeName],
trParam.thermo );
break;
case TP_ELECTCOND:
data.electCond = newLTP( xmlChild,
@ -1005,8 +1066,8 @@ namespace Cantera {
}
}
datatable[name] = data;
datatable.insert(pair<std::string, LiquidTransportData>(name,data));
// datatable[name] = data;
}
}
catch(CanteraError) {
@ -1059,12 +1120,10 @@ namespace Cantera {
for (int iChild = 0; iChild < num; iChild++) {
//tranTypeNode is a type of transport property like viscosity
XML_Node &tranTypeNode = transportNode.child(iChild);
if (tranTypeNode.hasChild("compositionDependence")) {
std::string nodeName = tranTypeNode.name();
if ( tranTypeNode.hasChild("compositionDependence")) {
//compDepNode contains the interaction model
XML_Node &compDepNode = tranTypeNode.child("compositionDependence");
std::string nodeName = tranTypeNode.name();
switch (m_tranPropMap[nodeName]) {
break;
case TP_VISCOSITY:
@ -1072,6 +1131,43 @@ namespace Cantera {
m_tranPropMap[nodeName],
trParam );
break;
case TP_IONCONDUCTIVITY:
trParam.ionConductivity = newLTI( compDepNode,
m_tranPropMap[nodeName],
trParam );
break;
case TP_MOBILITYRATIO:
{
int iSpec;
int numSpec = compDepNode.nChildren();
trParam.mobRatIndex.resize(numSpec);
trParam.mobilityRatio.resize(numSpec);
for (iSpec = 0; iSpec< numSpec; iSpec++){
XML_Node &propSpecNode = compDepNode.child(iSpec);
std::string specName = propSpecNode.name();
trParam.mobRatIndex[iSpec] = specName;
trParam.mobilityRatio[iSpec] = newLTI( propSpecNode,
m_tranPropMap[nodeName],
trParam );
};
};
break;
case TP_SELFDIFFUSION:
{
int iSpec;
int numSpec = compDepNode.nChildren();
trParam.selfDiffIndex.resize(numSpec);
trParam.selfDiffusion.resize(numSpec);
for (iSpec = 0; iSpec< numSpec; iSpec++){
XML_Node &propSpecNode = compDepNode.child(iSpec);
std::string specName = propSpecNode.name();
trParam.selfDiffIndex[iSpec] = specName;
trParam.selfDiffusion[iSpec] = newLTI( propSpecNode,
m_tranPropMap[nodeName],
trParam );
};
};
break;
case TP_THERMALCOND:
trParam.thermalCond = newLTI( compDepNode,
m_tranPropMap[nodeName],
@ -1122,15 +1218,17 @@ namespace Cantera {
trParam.velocityBasis_ = trParam.thermo->speciesIndex(velocityBasis) ;
else {
int linenum;
throw TransportDBError( linenum, "Unknown attribute " + velocityBasis + " for <velocityBasis> node. ");
throw TransportDBError( linenum, "Unknown attribute \"" + velocityBasis + "\" for <velocityBasis> node. ");
}
}
}
}
catch(CanteraError) {
;
catch (CanteraError) {
showErrors(std::cout);
}
//catch(CanteraError) {
// ;
//}
return;
}