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

*** empty log message ***

Browse source
This commit is contained in:
Dave Goodwin 2003-11-14 04:29:38 +00:00
parent 05b839c561
commit 2a64690ab7
4 changed files with 95 additions and 313 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

141
Cantera/src/transport/FtnTransport.h
View file

@ -7,14 +7,13 @@
* managers in Fortran.
*/
// Copyright 2001 California Institute of Technology
// Copyright 2003 California Institute of Technology
#ifndef CT_FTNTRANSPORT_H
#define CT_FTNTRANSPORT_H
#include "ct_defs.h"
#include "DenseMatrix.h"
#include "TransportBase.h"
/**
@ -31,116 +30,104 @@
#define __VISC__ visc_
#define __BULKVISC__ bvisc_
#define __TCON__ tcon_
#define __SPVISC__ spvisc_
#define __SPCOND__ spcond_
#define __SPFLUXES__ spfluxes_
#define __TDIFF__ tdiff_
#define __MULTIDIFF__ multidiff_
#define __MIXDIFF__ mixdiff_
#define __UPT__ updatet_
#define __UPC__ updatec_
#define __INIT__ init_
#define __SIGMA__ sigma_
#define __GETMOBILITIES__ getmobilities_
extern "C" {
doublereal __VISC__(doublereal* t, doublereal* p, doublereal* x);
doublereal __BULKVISC__(doublereal* t, doublereal* p, doublereal* x);
doublereal __TCON__(doublereal* t, doublereal* p, doublereal* x);
void __SPVISC__(doublereal* t, doublereal* p, doublereal* x, doublereal* visc);
void __SPCOND__(doublereal* t, doublereal* p, doublereal* x, doublereal* cond);
void __SPFLUXES__(doublereal* t, doublereal* p, doublereal* x,
integer* ndim, doublereal* gradt, integer* ldx, doublereal* gradx,
integer* ldf, doublereal* fluxes);
void __TDIFF__(doublereal* t, doublereal* p, doublereal* x, doublereal* dt);
void __MULTIDIFF__(doublereal* t, doublereal* p, doublereal* x,
integer* ld, doublereal* d);
void __MIXDIFF__(doublereal* t, doublereal* p, doublereal* x, doublereal* d);
void __BINDIFF__(doublereal* t, doublereal* p, doublereal* x, doublereal* d);
void __UPT__(doublereal* t);
void __UPC__(doublereal* p, doublereal* x);
void __INIT__();
void __BINDIFF__(doublereal* t, doublereal* p, doublereal* x, integer* ld, doublereal* d);
doublereal __SIGMA__(doublereal* t, doublereal* p, doublereal* x);
doublereal __GETMOBILITIES__(doublereal* t, doublereal* p,
doublereal* x, doublereal* mobil);
}
namespace Cantera {
/**
* A class that calls external Fortran functions to evaluate
* transport properties. Not currently used - may need updating.
* transport properties.
*/
class FtnTransport : public Transport {
public:
FtnTransport(int model) { m_model = model; }
FtnTransport(int model, thermo_t* thermo) : Transport(thermo) {
m_model = model;
m_x.resize(m_thermo->nSpecies(), 0.0);
updateTPX();
}
virtual int model() { return m_model; }
virtual int model() { return cFtnTransport + m_model; }
virtual doublereal viscosity() {
updateTPX();
return __VISC__(&m_temp, &m_pres, m_x.begin());
}
virtual void getSpeciesViscosities(doublereal* visc) {
__SPVISC__(&m_temp, &m_pres, m_x.begin(), visc);
}
virtual void getSpeciesConductivities(doublereal* cond) {
__SPCOND__(&m_temp, &m_pres, m_x.begin(), cond);
}
virtual doublereal bulkViscosity()
{ return __BULKVISC__(&m_temp, &m_pres, m_x.begin()); }
virtual doublereal thermalConductivity()
{ return __TCON__(&m_temp, &m_pres, m_x.begin()); }
virtual void getSpeciesFluxes(doublereal p, int ndim,
doublereal* grad_T, int ldx, doublereal* grad_X,
int ldf, doublereal* fluxes) {
doublereal pp = p;
integer ldxx = ldx, ndimm = ndim, ldff = ldf;
__SPFLUXES__(&m_temp, &pp, &m_x, &ndimm, grad_T, &ldxx, grad_X,
&ldff, fluxes);
virtual doublereal bulkViscosity() {
updateTPX();
return __BULKVISC__(&m_temp, &m_pres, m_x.begin());
}
virtual void getThermalDiffCoeffs(doublereal* dt)
{ __TDIFF__(&m_temp, &m_pres, m_x.begin(), dt); }
virtual doublereal thermalConductivity() {
updateTPX();
return __TCON__(&m_temp, &m_pres, m_x.begin());
}
virtual doublereal electricalConductivity() {
updateTPX();
return __SIGMA__(&m_temp, &m_pres, m_x.begin());
}
virtual void getMobilities(doublereal* mobil) {
updateTPX();
__GETMOBILITIES__(&m_temp, &m_pres, m_x.begin(), mobil);
}
virtual void getBinaryDiffCoeffs(doublereal p, int ld, doublereal* d)
{ m_pres = p;
integer ldd = ld;
__BINDIFF__(&m_temp, &m_pres, m_x.begin(), &ldd, d);
}
virtual void getMultiDiffCoeffs(doublereal p, int ld, doublereal* d)
{ m_pres = p;
integer ldd = ld;
__MULTIDIFF__(&m_temp, &m_pres, m_x.begin(), &ldd, d);
}
virtual void getMixDiffCoeffs(doublereal p, doublereal* d)
{ m_pres = p;
integer ldd = ld;
__MIXDIFF__(&m_temp, &m_pres, m_x.begin(), d);
}
virtual void update_T()
{ m_temp = m_mix->temperature();
__UPT__(m_temp);
}
virtual void update_C()
{ m_mix->getMoleFractions(m_x.begin());
m_pres = m_mix->pressure();
__UPC__(m_pres, m_x);
}
virtual bool init(TransportParams& tr) {
m_mix = tr.mix;
__INIT__();
virtual void getThermalDiffCoeffs(doublereal* dt) {
updateTPX();
__TDIFF__(&m_temp, &m_pres, m_x.begin(), dt);
}
virtual void getBinaryDiffCoeffs(int ld, doublereal* d) {
updateTPX();
integer ldd = ld;
__BINDIFF__(&m_temp, &m_pres, m_x.begin(), &ldd, d);
}
virtual void getMultiDiffCoeffs(int ld, doublereal* d) {
updateTPX();
integer ldd = ld;
__MULTIDIFF__(&m_temp, &m_pres, m_x.begin(), &ldd, d);
}
virtual void getMixDiffCoeffs(doublereal* d) {
updateTPX();
__MIXDIFF__(&m_temp, &m_pres, m_x.begin(), d);
}
private:
void updateTPX() {
m_temp = m_thermo->temperature();
m_pres = m_thermo->pressure();
m_thermo->getMoleFractions(m_x.begin());
}
doublereal m_temp;
doublereal m_pres;
vector_fp m_x;

9
Cantera/src/transport/MultiTransport.h
View file

@ -24,15 +24,6 @@
#endif
// STL includes
#include <vector>
#include <string>
#include <map>
#include <numeric>
#include <algorithm>
using namespace std;
// Cantera includes
#include "TransportBase.h"
#include "../DenseMatrix.h"

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

@ -49,6 +49,8 @@ namespace Cantera {
const int CK_MixtureAveraged = 211;
const int cSolidTransport = 300;
const int cDustyGasTransport = 400;
const int cUserTransport = 500;
const int cFtnTransport = 600;
class XML_Writer;

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

@ -19,6 +19,8 @@
#include "MixTransport.h"
#include "SolidTransport.h"
#include "DustyGasTransport.h"
#include "FtnTransport.h"
#include "TransportFactory.h"
#include "polyfit.h"
@ -30,8 +32,6 @@
#include "IdealGasPhase.h"
#include "ctml.h"
#include <iostream>
/**
* polynomial degree used for fitting collision integrals
@ -82,7 +82,7 @@ namespace Cantera {
* first-order binary diffusion coefficient should be multiplied
* to produce the value correct to second order. The expressions
* here are taken from Marerro and Mason,
* J. Phys. Chem. Ref. Data, vol. 1, p. 3 (1972).
* J. Phys. Chem. Ref. Data, vol. 1, p. 3 (1972).
*
* @param t Temperature (K)
* @param tr Transport parameters
@ -92,6 +92,8 @@ namespace Cantera {
* @param xmj mole fraction of species j
* @param fkj multiplier for d(k,j)
* @param fjk multiplier for d(j,k)
*
* @note This method is not used currently.
*/
void TransportFactory::getBinDiffCorrection(doublereal t,
const TransportParams& tr, int k, int j, doublereal xk, doublereal xj,
@ -168,6 +170,7 @@ namespace Cantera {
void TransportFactory::makePolarCorrections(int i, int j,
const TransportParams& tr, doublereal& f_eps, doublereal& f_sigma) {
// no correction if both are nonpolar, or both are polar
if (tr.polar[i] == tr.polar[j]) {
f_eps = 1.0; f_sigma = 1.0; return;
}
@ -182,8 +185,7 @@ namespace Cantera {
d3np = pow(tr.sigma[knp],3);
d3p = pow(tr.sigma[kp],3);
alpha_star = tr.alpha[knp]/d3np;
mu_p_star = tr.dipole(kp,kp)/sqrt(d3p
* tr.eps[kp]);
mu_p_star = tr.dipole(kp,kp)/sqrt(d3p * tr.eps[kp]);
xi = 1.0 + 0.25 * alpha_star * mu_p_star * mu_p_star *
sqrt(tr.eps[kp]/tr.eps[knp]);
f_sigma = pow(xi, -1.0/6.0);
@ -191,7 +193,6 @@ namespace Cantera {
}
// Constructor does nothing
TransportFactory::TransportFactory() : m_integrals(0) {
m_models["Mix"] = cMixtureAveraged;
m_models["Multi"] = cMulticomponent;
@ -199,6 +200,7 @@ namespace Cantera {
m_models["DustyGas"] = cDustyGasTransport;
m_models["CK_Multi"] = CK_Multicomponent;
m_models["CK_Mix"] = CK_MixtureAveraged;
m_models["User"] = cUserTransport;
m_models["None"] = 0;
}
@ -206,22 +208,21 @@ namespace Cantera {
/**
* Destructor
*
* We do not delete statically
* created single instance of this class here, because it would
* create an infinite loop if destructor is called for that
* single instance.
* However, we do have a malloced pointer to m_integrals
* that does need to be explicitly deleted.
* We do not delete statically created single instance of this
* class here, because it would create an infinite loop if
* destructor is called for that single instance. However, we do
* have a pointer to m_integrals that does need to be
* explicitly deleted.
*/
TransportFactory::~TransportFactory() {
if (m_integrals) {
delete m_integrals;
m_integrals = 0;
delete m_integrals;
m_integrals = 0;
}
}
/**
* This static function deletes the statically malloced instance.
* This static function deletes the statically allocated instance.
*/
void TransportFactory::deleteTransportFactory() {
if (__factory) {
@ -281,35 +282,6 @@ namespace Cantera {
return tr;
}
// Transport* TransportFactory::
// newPowerTransport(const string& transport_database,
// Transport::phase_t* mix) {
// doublereal tref, viscref, lambdaref, viscexp, lambdaexp,
// diffexp;
// ifstream f(transport_database.c_str());
// if (!f) {
// throw TransportDBError("newPowerTransport: file" +
// transport_database + " not found.");
// }
// f >> tref
// >> viscref >> viscexp
// >> lambdaref >> lambdaexp
// >> diffexp;
// int nsp = mix->nSpecies();
// int k;
// vector_fp dref(nsp);
// for (k = 0; k < nsp; k++) {
// f >> dref[k];
// }
// PowerTransport* t = new PowerTransport();
// t->init(mix);
// t->setExponents(viscexp, lambdaexp, diffexp);
// t->setRefValues(tref, viscref, lambdaref, nsp, dref.begin());
// f.close();
// return t;
// }
/**
* Prepare to build a new kinetic-theory-based transport manager
@ -321,7 +293,6 @@ namespace Cantera {
thermo_t* thermo, int mode, int log_level, TransportParams& tr) {
// constant mixture attributes
//tr.mix = thermo;
tr.thermo = thermo;
tr.nsp = tr.thermo->nSpecies();
int nsp = tr.nsp;
@ -348,14 +319,7 @@ namespace Cantera {
tr.sigma.resize(nsp);
tr.eps.resize(nsp);
//readTransportDatabase(flog, transport_database,
// tr.thermo->speciesNames(), tr);
XML_Node root, log;
//string infile = findInputFile(transport_database);
//ifstream fin(infile.c_str());
//root.build(fin);
getTransportData(transport_database, log,
tr.thermo->speciesNames(), tr);
@ -413,10 +377,13 @@ namespace Cantera {
// Chemkin fits the entire T* range in the Monchick and Mason tables,
// so modify tstar_min and tstar_max if in Chemkin compatibility mode
//if (mode == CK_Mode) {
// NOTE: the 'if' was commented out DGG 11/12/03
if (mode == CK_Mode) { // uncommented
tstar_min = 0.101;
tstar_max = 99.9;
//}
} // uncommented
// initialize the collision integral calculator for the desired
// T* range
@ -540,152 +507,7 @@ namespace Cantera {
*
*********************************************************/
/**
* Read transport property data from a file for a list of species.
* Given the name of a file containing transport property
* parameters and a list of species names, this method returns an
* instance of TransportParams containing the transport data for
* these species read from the file.
*/
// void TransportFactory::readTransportDatabase(
// ostream& logfile,
// const string& database_file,
// const vector<string>& names,
// TransportParams& tr)
// {
// string dbase_file = findInputFile(database_file);
// ifstream dbase(dbase_file.c_str());
// if (!dbase) throw CanteraError("readTransportDatabase",
// "cannot open file " + database_file);
// string name, rest;
// int geom, linenum = 0;
// map<string, TransportData> datatable;
// doublereal welldepth, diam, dipole, polar, rot;
// // read all entries in database into 'datatable' and check for
// // errors. Note that this procedure validates all entries, not
// // only those for the species listed in 'names'.
// while (!dbase.eof()) {
// dbase >> name;
// linenum++;
// if (name[0] != '!' && !dbase.eof()) {
// dbase >> geom >> welldepth >> diam
// >> dipole >> polar >> rot;
// // get the rest of the line, in case there are comments
// getline(dbase, rest);
// TransportData data;
// data.speciesName = name;
// data.geometry = geom;
// data.wellDepth = welldepth;
// data.diameter = diam;
// data.dipoleMoment = dipole;
// data.polarizability = polar;
// data.rotRelaxNumber = rot;
// datatable[name] = data;
// if (welldepth >= 0.0) data.wellDepth = welldepth;
// else throw TransportDBError(linenum,
// "negative well depth");
// if (diam > 0.0) data.diameter = diam;
// else throw TransportDBError(linenum,
// "negative or zero diameter");
// if (dipole >= 0.0) data.dipoleMoment = dipole;
// else throw TransportDBError(linenum,
// "negative dipole moment");
// if (polar >= 0.0) data.polarizability = polar;
// else throw TransportDBError(linenum,
// "negative polarizability");
// if (rot >= 0.0) data.rotRelaxNumber = rot;
// else throw TransportDBError(linenum,
// "negative rotation relaxation number");
// datatable[name] = data;
// }
// }
// dbase.close();
// // look up the entries for the species listed in 'names'.
// tr.xml->XML_open(logfile,"database_parameters");
// tr.xml->XML_item(logfile,"file_name",database_file);
// for (int i = 0; i < tr.nsp; i++) {
// TransportData& trdat = datatable[names[i]];
// // 'datatable' returns a default TransportData object if
// // the species name is not one in the transport database.
// // This can be detected by examining 'geometry'.
// if (trdat.geometry < 0) {
// throw TransportDBError(0,"no transport data found for species "
// + names[i]);
// }
// // parameters are converted to SI units before storing
// // rotational heat capacity / R
// switch (trdat.geometry) {
// case 0:
// tr.crot[i] = 0.0; // monatomic
// break;
// case 1:
// tr.crot[i] = 1.0; // linear
// break;
// default:
// tr.crot[i] = 1.5; // nonlinear
// }
// tr.dipole(i,i) = 1.e-25 * SqrtTen * trdat.dipoleMoment;
// if (trdat.dipoleMoment > 0.0)
// tr.polar[i] = true;
// else
// tr.polar[i] = false;
// // A^3 -> m^3
// tr.alpha[i] = 1.e-30 * trdat.polarizability;
// tr.sigma[i] = 1.e-10 * trdat.diameter;
// tr.eps[i] = Boltzmann * trdat.wellDepth;
// tr.zrot[i] = fmaxx(1.0, trdat.rotRelaxNumber);
// // write database parameters to log file
// tr.xml->XML_open(logfile, names[i]);
// tr.xml->XML_item(logfile, "geom", trdat.geometry);
// tr.xml->XML_item(logfile, "epsilon", trdat.wellDepth);
// tr.xml->XML_item(logfile, "sigma", trdat.diameter);
// tr.xml->XML_item(logfile, "dipole", trdat.dipoleMoment);
// tr.xml->XML_item(logfile, "alpha", trdat.polarizability);
// tr.xml->XML_item(logfile, "zrot", trdat.rotRelaxNumber);
// tr.xml->XML_close(logfile, names[i]);
// }
// tr.xml->XML_close(logfile,"database_parameters");
// }
/*********************************************************
*
* Read Transport Database
*
*********************************************************/
/**
/**
* Read transport property data from a file for a list of species.
* Given the name of a file containing transport property
* parameters and a list of species names, this method returns an
@ -730,14 +552,7 @@ namespace Cantera {
polar = getFloat(tr, "polarizability");
rot = getFloat(tr, "rotRelax");
//getFloats(tr, fv, false);
//welldepth = fv["LJ_welldepth"];
//diam = fv["LJ_diameter"];
//dipole = fv["dipoleMoment"];
//polar = fv["polarizability"];
//rot = fv["rotRelax"];
GasTransportData data;
GasTransportData data;
data.speciesName = name;
data.geometry = geom;
if (welldepth >= 0.0) data.wellDepth = welldepth;
@ -763,11 +578,7 @@ namespace Cantera {
datatable[name] = data;
}
// look up the entries for the species listed in 'names'.
//tr.xml->XML_open(logfile,"database_parameters");
//tr.xml->XML_item(logfile,"file_name",database_file);
for (i = 0; i < tr.nsp; i++) {
for (i = 0; i < tr.nsp; i++) {
GasTransportData& trdat = datatable[names[i]];
@ -809,18 +620,8 @@ namespace Cantera {
tr.eps[i] = Boltzmann * trdat.wellDepth;
tr.zrot[i] = fmaxx(1.0, trdat.rotRelaxNumber);
// write database parameters to log file
//tr.xml->XML_open(logfile, names[i]);
//tr.xml->XML_item(logfile, "geom", trdat.geometry);
//tr.xml->XML_item(logfile, "epsilon", trdat.wellDepth);
//tr.xml->XML_item(logfile, "sigma", trdat.diameter);
//tr.xml->XML_item(logfile, "dipole", trdat.dipoleMoment);
//tr.xml->XML_item(logfile, "alpha", trdat.polarizability);
//tr.xml->XML_item(logfile, "zrot", trdat.rotRelaxNumber);
//tr.xml->XML_close(logfile, names[i]);
}
//tr.xml->XML_close(logfile,"database_parameters");
}
}
}
/*********************************************************
@ -902,7 +703,6 @@ namespace Cantera {
tr.xml->XML_comment(logfile,s);
}
// const vector_fp& cp_R = tr.mix->cp_R();
doublereal cp_R, cond, w_RT, f_int, A_factor, B_factor,
c1, cv_rot, cv_int, f_rot, f_trans, om11;
@ -1072,13 +872,15 @@ namespace Cantera {
pow((Boltzmann * t), 1.5)/
(Pi * sigma * sigma * om11);
// 2nd order correction
// NOTE: THIS CORRECTION IS NOT APPLIED
doublereal fkj, fjk;
getBinDiffCorrection(t, tr, k, j, 1.0, 1.0, fkj, fjk);
//diffcoeff *= fkj;
if (mode == CK_Mode) {
if (mode == CK_Mode) {
diff[n] = log(diffcoeff);
w[n] = -1.0;
}