Doxygen updates, and code formatting. No actual changes.
This commit is contained in:
parent
0bbd0879d5
commit
900f5d1530
2 changed files with 226 additions and 111 deletions
|
|
@ -761,47 +761,45 @@ namespace Cantera {
|
|||
tr.xml->XML_comment(logfile, "*** polynomial coefficients not printed (log_level < 3) ***");
|
||||
}
|
||||
#endif
|
||||
for (i = 0; i < nsp; i++)
|
||||
{
|
||||
for (j = i; j < nsp; j++)
|
||||
{
|
||||
// Chemkin fits only delta* = 0
|
||||
if (mode != CK_Mode)
|
||||
dstar = tr.delta(i,j);
|
||||
else
|
||||
dstar = 0.0;
|
||||
for (i = 0; i < nsp; i++) {
|
||||
for (j = i; j < nsp; j++) {
|
||||
// Chemkin fits only delta* = 0
|
||||
if (mode != CK_Mode) {
|
||||
dstar = tr.delta(i,j);
|
||||
} else {
|
||||
dstar = 0.0;
|
||||
}
|
||||
|
||||
// if a fit has already been generated for
|
||||
// delta* = tr.delta(i,j), then use it. Otherwise,
|
||||
// make a new fit, and add tr.delta(i,j) to the list
|
||||
// of delta* values for which fits have been done.
|
||||
// if a fit has already been generated for
|
||||
// delta* = tr.delta(i,j), then use it. Otherwise,
|
||||
// make a new fit, and add tr.delta(i,j) to the list
|
||||
// of delta* values for which fits have been done.
|
||||
|
||||
// 'find' returns a pointer to end() if not found
|
||||
if (dptr = find(tr.fitlist.begin(), tr.fitlist.end(),
|
||||
dstar), dptr == tr.fitlist.end())
|
||||
{
|
||||
vector_fp ca(degree+1), cb(degree+1), cc(degree+1);
|
||||
vector_fp co22(degree+1);
|
||||
m_integrals->fit(logfile, degree, dstar,
|
||||
DATA_PTR(ca), DATA_PTR(cb), DATA_PTR(cc));
|
||||
m_integrals->fit_omega22(logfile, degree, dstar,
|
||||
DATA_PTR(co22));
|
||||
tr.omega22_poly.push_back(co22);
|
||||
tr.astar_poly.push_back(ca);
|
||||
tr.bstar_poly.push_back(cb);
|
||||
tr.cstar_poly.push_back(cc);
|
||||
tr.poly[i][j] = static_cast<int>(tr.astar_poly.size()) - 1;
|
||||
tr.fitlist.push_back(dstar);
|
||||
}
|
||||
// 'find' returns a pointer to end() if not found
|
||||
dptr = find(tr.fitlist.begin(), tr.fitlist.end(), dstar);
|
||||
if (dptr == tr.fitlist.end()) {
|
||||
vector_fp ca(degree+1), cb(degree+1), cc(degree+1);
|
||||
vector_fp co22(degree+1);
|
||||
m_integrals->fit(logfile, degree, dstar,
|
||||
DATA_PTR(ca), DATA_PTR(cb), DATA_PTR(cc));
|
||||
m_integrals->fit_omega22(logfile, degree, dstar,
|
||||
DATA_PTR(co22));
|
||||
tr.omega22_poly.push_back(co22);
|
||||
tr.astar_poly.push_back(ca);
|
||||
tr.bstar_poly.push_back(cb);
|
||||
tr.cstar_poly.push_back(cc);
|
||||
tr.poly[i][j] = static_cast<int>(tr.astar_poly.size()) - 1;
|
||||
tr.fitlist.push_back(dstar);
|
||||
}
|
||||
|
||||
// delta* found in fitlist, so just point to this
|
||||
// polynomial
|
||||
else {
|
||||
tr.poly[i][j] = static_cast<int>((dptr - tr.fitlist.begin()));
|
||||
}
|
||||
tr.poly[j][i] = tr.poly[i][j];
|
||||
}
|
||||
}
|
||||
// delta* found in fitlist, so just point to this
|
||||
// polynomial
|
||||
else {
|
||||
tr.poly[i][j] = static_cast<int>((dptr - tr.fitlist.begin()));
|
||||
}
|
||||
tr.poly[j][i] = tr.poly[i][j];
|
||||
}
|
||||
}
|
||||
#ifdef DEBUG_MODE
|
||||
if (m_verbose) {
|
||||
tr.xml->XML_close(logfile, "tstar_fits");
|
||||
|
|
@ -1462,71 +1460,69 @@ namespace Cantera {
|
|||
mxerr = 0.0, mxrelerr = 0.0;
|
||||
vector_fp diff(np + 1);
|
||||
doublereal eps, sigma;
|
||||
for (k = 0; k < tr.nsp_; k++)
|
||||
{
|
||||
for (j = k; j < tr.nsp_; j++) {
|
||||
for (k = 0; k < tr.nsp_; k++) {
|
||||
for (j = k; j < tr.nsp_; j++) {
|
||||
|
||||
ipoly = tr.poly[k][j];
|
||||
for (n = 0; n < np; n++) {
|
||||
ipoly = tr.poly[k][j];
|
||||
for (n = 0; n < np; n++) {
|
||||
|
||||
t = tr.tmin + dt*n;
|
||||
t = tr.tmin + dt*n;
|
||||
|
||||
eps = tr.epsilon(j,k);
|
||||
tstar = Boltzmann * t/eps;
|
||||
sigma = tr.diam(j,k);
|
||||
om11 = m_integrals->omega11(tstar, tr.delta(j,k));
|
||||
eps = tr.epsilon(j,k);
|
||||
tstar = Boltzmann * t/eps;
|
||||
sigma = tr.diam(j,k);
|
||||
om11 = m_integrals->omega11(tstar, tr.delta(j,k));
|
||||
|
||||
diffcoeff = ThreeSixteenths *
|
||||
sqrt( 2.0 * Pi/tr.reducedMass(k,j) ) *
|
||||
pow((Boltzmann * t), 1.5)/
|
||||
(Pi * sigma * sigma * om11);
|
||||
diffcoeff = ThreeSixteenths *
|
||||
sqrt( 2.0 * Pi/tr.reducedMass(k,j) ) *
|
||||
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;
|
||||
// 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) {
|
||||
diff[n] = log(diffcoeff);
|
||||
w[n] = -1.0;
|
||||
}
|
||||
else {
|
||||
diff[n] = diffcoeff/pow(t, 1.5);
|
||||
w[n] = 1.0/(diff[n]*diff[n]);
|
||||
}
|
||||
if (mode == CK_Mode) {
|
||||
diff[n] = log(diffcoeff);
|
||||
w[n] = -1.0;
|
||||
}
|
||||
polyfit(np, DATA_PTR(tlog), DATA_PTR(diff),
|
||||
DATA_PTR(w), degree, ndeg, 0.0, DATA_PTR(c));
|
||||
|
||||
doublereal pre;
|
||||
for (n = 0; n < np; n++) {
|
||||
if (mode == CK_Mode) {
|
||||
val = exp(diff[n]);
|
||||
fit = exp(poly3(tlog[n], DATA_PTR(c)));
|
||||
}
|
||||
else {
|
||||
t = exp(tlog[n]);
|
||||
pre = pow(t, 1.5);
|
||||
val = pre * diff[n];
|
||||
fit = pre * poly4(tlog[n], DATA_PTR(c));
|
||||
}
|
||||
err = fit - val;
|
||||
relerr = err/val;
|
||||
if (fabs(err) > mxerr) mxerr = fabs(err);
|
||||
if (fabs(relerr) > mxrelerr) mxrelerr = fabs(relerr);
|
||||
else {
|
||||
diff[n] = diffcoeff/pow(t, 1.5);
|
||||
w[n] = 1.0/(diff[n]*diff[n]);
|
||||
}
|
||||
tr.diffcoeffs.push_back(c);
|
||||
#ifdef DEBUG_MODE
|
||||
if (tr.log_level >= 2 && m_verbose) {
|
||||
tr.xml->XML_writeVector(logfile, " ", tr.thermo->speciesName(k)
|
||||
+ "__"+tr.thermo->speciesName(j), c.size(), DATA_PTR(c));
|
||||
}
|
||||
#endif
|
||||
}
|
||||
polyfit(np, DATA_PTR(tlog), DATA_PTR(diff),
|
||||
DATA_PTR(w), degree, ndeg, 0.0, DATA_PTR(c));
|
||||
|
||||
doublereal pre;
|
||||
for (n = 0; n < np; n++) {
|
||||
if (mode == CK_Mode) {
|
||||
val = exp(diff[n]);
|
||||
fit = exp(poly3(tlog[n], DATA_PTR(c)));
|
||||
} else {
|
||||
t = exp(tlog[n]);
|
||||
pre = pow(t, 1.5);
|
||||
val = pre * diff[n];
|
||||
fit = pre * poly4(tlog[n], DATA_PTR(c));
|
||||
}
|
||||
err = fit - val;
|
||||
relerr = err/val;
|
||||
if (fabs(err) > mxerr) mxerr = fabs(err);
|
||||
if (fabs(relerr) > mxrelerr) mxrelerr = fabs(relerr);
|
||||
}
|
||||
tr.diffcoeffs.push_back(c);
|
||||
#ifdef DEBUG_MODE
|
||||
if (tr.log_level >= 2 && m_verbose) {
|
||||
tr.xml->XML_writeVector(logfile, " ", tr.thermo->speciesName(k)
|
||||
+ "__"+tr.thermo->speciesName(j), c.size(), DATA_PTR(c));
|
||||
}
|
||||
#endif
|
||||
}
|
||||
}
|
||||
#ifdef DEBUG_MODE
|
||||
if (m_verbose) {
|
||||
sprintf(s,"Maximum binary diffusion coefficient absolute error:"
|
||||
|
|
|
|||
|
|
@ -86,12 +86,12 @@ namespace Cantera {
|
|||
|
||||
public:
|
||||
|
||||
//! Constructor
|
||||
GasTransportParams() :
|
||||
TransportParams(),
|
||||
visccoeffs(0),
|
||||
condcoeffs(0),
|
||||
diffcoeffs(0),
|
||||
polytempvec(0),
|
||||
poly(0),
|
||||
omega22_poly(0),
|
||||
astar_poly(0),
|
||||
|
|
@ -112,32 +112,89 @@ namespace Cantera {
|
|||
{
|
||||
}
|
||||
|
||||
//! Destructor
|
||||
virtual ~GasTransportParams() {}
|
||||
|
||||
// polynomial fits
|
||||
|
||||
//! temperature-fit of the viscosity
|
||||
/*!
|
||||
* The outer loop the number of species, nsp
|
||||
* The inner loop is over degree + 1, which is the polynomial order of the collision integral fit.
|
||||
*/
|
||||
std::vector<vector_fp> visccoeffs;
|
||||
|
||||
//temperature-fit heat conduction
|
||||
std::vector<vector_fp> condcoeffs;
|
||||
//! temperature-fits of the heat conduction
|
||||
/*!
|
||||
* The outer loop the number of species, nsp
|
||||
* The inner loop is over degree + 1, which is the polynomial order of the collision integral fit.
|
||||
*/
|
||||
std::vector<vector_fp> condcoeffs;
|
||||
|
||||
//temperature-fit diffusivity
|
||||
std::vector<vector_fp> diffcoeffs;
|
||||
vector_fp polytempvec;
|
||||
//! temperature-fits of the diffusivity
|
||||
/*!
|
||||
* The outer loop the number of species, nsp
|
||||
* The inner loop is over degree + 1, which is the polynomial order of the collision integral fit.
|
||||
*/
|
||||
std::vector<vector_fp> diffcoeffs;
|
||||
|
||||
|
||||
//! This is vector of vectors containing the integer ookup value for the (i,j) interaction
|
||||
/*!
|
||||
* The outer loop is over a flat (i,j) index that is parameterized on the tr.delta(i,j) value.
|
||||
* Unique values of delta get their own spot in the array. The values of delta are storred in
|
||||
* the fitlist vector.
|
||||
*
|
||||
* The inner loop is over degree + 1, which is the polynomial order of the collision integral fit.
|
||||
*/
|
||||
std::vector<std::vector<int> > poly;
|
||||
std::vector<vector_fp > omega22_poly;
|
||||
std::vector<vector_fp > astar_poly;
|
||||
std::vector<vector_fp > bstar_poly;
|
||||
std::vector<vector_fp > cstar_poly;
|
||||
|
||||
//! This is vector of vectors containing the astar fit.
|
||||
/*!
|
||||
* The outer loop is over a flat (i,j) index that is parameterized on the tr.delta(i,j) value.
|
||||
* Unique values of delta get their own spot in the array. The values of delta are storred in
|
||||
* the fitlist vector.
|
||||
*
|
||||
* The inner loop is over degree + 1, which is the polynomial order of the collision integral fit.
|
||||
*/
|
||||
std::vector<vector_fp> omega22_poly;
|
||||
|
||||
//! This is vector of vectors containing the astar fit.
|
||||
/*!
|
||||
* The outer loop is over a flat (i,j) index that is parameterized on the tr.delta(i,j) value.
|
||||
* Unique values of delta get their own spot in the array. The values of delta are storred in
|
||||
* the fitlist vector.
|
||||
*
|
||||
* The inner loop is over degree + 1, which is the polynomial order of the collision integral fit.
|
||||
*/
|
||||
std::vector<vector_fp> astar_poly;
|
||||
|
||||
//! This is vector of vectors containing the astar fit.
|
||||
/*!
|
||||
* The outer loop is over a flat (i,j) index that is parameterized on the tr.delta(i,j) value.
|
||||
* Unique values of delta get their own spot in the array. The values of delta are storred in
|
||||
* the fitlist vector.
|
||||
*
|
||||
* The inner loop is over degree + 1, which is the polynomial order of the collision integral fit.
|
||||
*/
|
||||
std::vector<vector_fp> bstar_poly;
|
||||
|
||||
//! This is vector of vectors containing the astar fit.
|
||||
/*!
|
||||
* The outer loop is over a flat (i,j) index that is parameterized on the tr.delta(i,j) value.
|
||||
* Unique values of delta get their own spot in the array. The values of delta are storred in
|
||||
* the fitlist vector.
|
||||
*
|
||||
* The inner loop is over degree + 1, which is the polynomial order of the collision integral fit.
|
||||
*/
|
||||
std::vector<vector_fp> cstar_poly;
|
||||
|
||||
//! Rotational relaxation number for the species in the current phase
|
||||
/*!
|
||||
* length is the number of species in the phase
|
||||
* units are dimensionless
|
||||
*/
|
||||
vector_fp zrot;
|
||||
vector_fp zrot;
|
||||
|
||||
//! Dimensionless rotational heat capacity of the species in the current phase
|
||||
/*!
|
||||
|
|
@ -145,31 +202,93 @@ namespace Cantera {
|
|||
* length is the number of species in the pahse
|
||||
* units are dimensionless (Cr / R)
|
||||
*/
|
||||
vector_fp crot;
|
||||
vector_fp crot;
|
||||
|
||||
//! Vector of booleans indicating whether a species is a polar molecule
|
||||
/*!
|
||||
* Length is nsp
|
||||
*/
|
||||
std::vector<bool> polar;
|
||||
|
||||
//! Polarizability of each species in the phase
|
||||
/*!
|
||||
* Length = nsp
|
||||
* Units = m^3
|
||||
*/
|
||||
vector_fp alpha;
|
||||
vector_fp fitlist;
|
||||
|
||||
//! This is vector containing the values of delta(i,j) that are used in the collision integral fits.
|
||||
/*!
|
||||
* This is used in astar_poly, bstar_poly, cstar_poly, and omega22_poly.
|
||||
* The outer loop is over a flat (i,j) index that is parameterized on the tr.delta(i,j) value.
|
||||
* Unique values of delta get their own spot in the array. The values of delta are storred in
|
||||
* the fitlist vector.
|
||||
*
|
||||
*/
|
||||
vector_fp fitlist;
|
||||
|
||||
//! Lennard-Jones well-depth of the species in the current phase
|
||||
/*!
|
||||
* length is the number of species in the phase
|
||||
* Units are Joules (Note this is not Joules/kmol)
|
||||
* Units are Joules (Note this is not Joules/kmol) (note, no kmol -> this is a per molecule amount)
|
||||
*/
|
||||
vector_fp eps;
|
||||
vector_fp eps;
|
||||
|
||||
//! Lennard-Jones diameter of the species in the current phase
|
||||
/*!
|
||||
* length is the number of species in the phase
|
||||
* units are in meters.
|
||||
*/
|
||||
vector_fp sigma;
|
||||
vector_fp sigma;
|
||||
|
||||
//! This is the reduced mass of the interaction between species i and j
|
||||
/*!
|
||||
* tr.reducedMass(i,j) = tr.mw[i] * tr.mw[j] / (Avogadro * (tr.mw[i] + tr.mw[j]));
|
||||
*
|
||||
* Units are kg (note, no kmol -> this is a per molecule amount)
|
||||
*
|
||||
* Length nsp * nsp. This is a symmetric matrix
|
||||
*/
|
||||
DenseMatrix reducedMass;
|
||||
DenseMatrix diam;
|
||||
DenseMatrix epsilon;
|
||||
DenseMatrix dipole;
|
||||
DenseMatrix delta;
|
||||
|
||||
//! hard-sphere diameter for (i,j) collision
|
||||
/*!
|
||||
* diam(i,j) = 0.5*(tr.sigma[i] + tr.sigma[j]);
|
||||
* Units are m (note, no kmol -> this is a per molecule amount)
|
||||
*
|
||||
* Length nsp * nsp. This is a symmetric matrix.
|
||||
*/
|
||||
DenseMatrix diam;
|
||||
|
||||
//! The effective well depth for (i,j) collisions
|
||||
/*!
|
||||
* epsilon(i,j) = sqrt(tr.eps[i]*tr.eps[j]);
|
||||
* Units are Joules (note, no kmol -> this is a per molecule amount)
|
||||
*
|
||||
* Length nsp * nsp. This is a symmetric matrix.
|
||||
*/
|
||||
DenseMatrix epsilon;
|
||||
|
||||
//! The effective dipole moment for (i,j) collisions
|
||||
/*!
|
||||
* tr.dipoleMoment has units of Debye's. A Debye is 10-18 cm3/2 erg1/2
|
||||
*
|
||||
* tr.dipole(i,i) = 1.e-25 * SqrtTen * trdat.dipoleMoment;
|
||||
* tr.dipole(i,j) = sqrt(tr.dipole(i,i)*tr.dipole(j,j));
|
||||
* Units are in Debye (note, no kmol -> this is a per molecule amount)
|
||||
*
|
||||
* Length nsp * nsp. This is a symmetric matrix.
|
||||
*/
|
||||
DenseMatrix dipole;
|
||||
|
||||
//! Matrix containing the reduced dipole moment of the interaction between two species
|
||||
/*!
|
||||
* This is the reduced dipole moment of the interaction between two species
|
||||
* 0.5 * tr.dipole(i,j)*tr.dipole(i,j) (epsilon(i,j) * d * d * d);
|
||||
*
|
||||
* Length nsp * nsp .This is a symmetric matrix
|
||||
*/
|
||||
DenseMatrix delta;
|
||||
|
||||
};
|
||||
|
||||
|
|
|
|||
Loading…
Add table
Reference in a new issue