[Transport] Move implementation of MultiTransport into a single file

This commit is contained in:
Ray Speth 2014-01-06 18:45:57 +00:00
parent 907bfa0a10
commit fc008d067d
2 changed files with 196 additions and 213 deletions

View file

@ -1,212 +0,0 @@
/**
* @file L_matrix.h
* Functions to evaluate portions of the L matrix needed for
* multicomponent transport properties.
*/
#ifndef CT_LMATRIX_H
#define CT_LMATRIX_H
#include "cantera/transport/MultiTransport.h"
namespace Cantera
{
//! Constant to compare dimensionless heat capacities against zero
const doublereal Min_C_Internal = 0.001;
bool MultiTransport::hasInternalModes(size_t j)
{
return (m_cinternal[j] > Min_C_Internal);
}
void MultiTransport::eval_L0000(const doublereal* const x)
{
doublereal prefactor = 16.0*m_temp/25.0;
doublereal sum;
for (size_t i = 0; i < m_nsp; i++) {
// subtract-off the k=i term to account for the first delta
// function in Eq. (12.121)
sum = -x[i]/m_bdiff(i,i);
for (size_t k = 0; k < m_nsp; k++) {
sum += x[k]/m_bdiff(i,k);
}
sum /= m_mw[i];
for (size_t j = 0; j != m_nsp; ++j) {
m_Lmatrix(i,j) = prefactor * x[j]
* (m_mw[j] * sum + x[i]/m_bdiff(i,j));
}
// diagonal term is zero
m_Lmatrix(i,i) = 0.0;
}
}
void MultiTransport::eval_L0010(const doublereal* const x)
{
doublereal prefactor = 1.6*m_temp;
doublereal sum, wj, xj;
for (size_t j = 0; j < m_nsp; j++) {
//constant = prefactor * x[j];
xj = x[j];
wj = m_mw[j];
sum = 0.0;
for (size_t i = 0; i < m_nsp; i++) {
m_Lmatrix(i,j + m_nsp) = - prefactor * x[i] * xj * m_mw[i] *
(1.2 * m_cstar(j,i) - 1.0) /
((wj + m_mw[i]) * m_bdiff(j,i));
// the next term is independent of "j";
// need to do it for the "j,j" term
sum -= m_Lmatrix(i,j+m_nsp);
}
m_Lmatrix(j,j+m_nsp) += sum;
}
}
void MultiTransport::eval_L1000()
{
for (size_t j = 0; j < m_nsp; j++) {
for (size_t i = 0; i < m_nsp; i++) {
m_Lmatrix(i+m_nsp,j) = m_Lmatrix(j,i+m_nsp);
}
}
}
void MultiTransport::eval_L1010(const doublereal* x)
{
const doublereal fiveover3pi = 5.0/(3.0*Pi);
doublereal prefactor = (16.0*m_temp)/25.0;
doublereal constant1, wjsq, constant2, constant3, constant4,
fourmj, threemjsq, sum, sumwij;;
doublereal term1, term2;
for (size_t j = 0; j < m_nsp; j++) {
// get constant terms that depend on just species "j"
constant1 = prefactor*x[j];
wjsq = m_mw[j]*m_mw[j];
constant2 = 13.75*wjsq;
constant3 = m_crot[j]/m_rotrelax[j];
constant4 = 7.5*wjsq;
fourmj = 4.0*m_mw[j];
threemjsq = 3.0*m_mw[j]*m_mw[j];
sum = 0.0;
for (size_t i = 0; i < m_nsp; i++) {
sumwij = m_mw[i] + m_mw[j];
term1 = m_bdiff(i,j) * sumwij*sumwij;
term2 = fourmj*m_astar(i,j)*(1.0 + fiveover3pi*
(constant3 +
(m_crot[i]/m_rotrelax[i]))); // see Eq. (12.125)
m_Lmatrix(i+m_nsp,j+m_nsp) = constant1*x[i]*m_mw[i] /(m_mw[j]*term1) *
(constant2 - threemjsq*m_bstar(i,j)
- term2*m_mw[j]);
sum += x[i] /(term1) *
(constant4 + m_mw[i]*m_mw[i]*
(6.25 - 3.0*m_bstar(i,j)) + term2*m_mw[i]);
}
m_Lmatrix(j+m_nsp,j+m_nsp) -= sum*constant1;
}
}
void MultiTransport::eval_L1001(const doublereal* x)
{
doublereal prefactor = 32.00*m_temp/(5.00*Pi);
doublereal constant, sum;
size_t n2 = 2*m_nsp;
int npoly = 0;
for (size_t j = 0; j < m_nsp; j++) {
// collect terms that depend only on "j"
if (hasInternalModes(j)) {
constant = prefactor*m_mw[j]*x[j]*m_crot[j]/(m_cinternal[j]*m_rotrelax[j]);
sum = 0.0;
for (size_t i = 0; i < m_nsp; i++) {
// see Eq. (12.127)
m_Lmatrix(i+m_nsp,j+n2) = constant * m_astar(j,i) * x[i] /
((m_mw[j] + m_mw[i]) * m_bdiff(j,i));
sum += m_Lmatrix(i+m_nsp,j+n2);
}
npoly++;
m_Lmatrix(j+m_nsp,j+n2) += sum;
} else {
for (size_t i = 0; i < m_nsp; i++) {
m_Lmatrix(i+m_nsp,j+n2) = 0.0;
}
}
}
}
void MultiTransport::eval_L0001()
{
size_t n2 = 2*m_nsp;
for (size_t j = 0; j < m_nsp; j++) {
for (size_t i = 0; i < m_nsp; i++) {
m_Lmatrix(i,j+n2) = 0.0;
}
}
}
void MultiTransport::eval_L0100()
{
size_t n2 = 2*m_nsp;
for (size_t j = 0; j < m_nsp; j++)
for (size_t i = 0; i < m_nsp; i++) {
m_Lmatrix(i+n2,j) = 0.0; // see Eq. (12.123)
}
}
void MultiTransport::eval_L0110()
{
size_t n2 = 2*m_nsp;
for (size_t j = 0; j < m_nsp; j++)
for (size_t i = 0; i < m_nsp; i++) {
m_Lmatrix(i+n2,j+m_nsp) = m_Lmatrix(j+m_nsp,i+n2); // see Eq. (12.123)
}
}
void MultiTransport::eval_L0101(const doublereal* x)
{
const doublereal fivepi = 5.00*Pi;
const doublereal eightoverpi = 8.0 / Pi;
doublereal prefactor = 4.00*m_temp;
size_t n2 = 2*m_nsp;
doublereal constant1, constant2, diff_int, sum;
for (size_t i = 0; i < m_nsp; i++) {
if (hasInternalModes(i)) {
// collect terms that depend only on "i"
constant1 = prefactor*x[i]/m_cinternal[i];
constant2 = 12.00*m_mw[i]*m_crot[i] /
(fivepi*m_cinternal[i]*m_rotrelax[i]);
sum = 0.0;
for (size_t k = 0; k < m_nsp; k++) {
// see Eq. (12.131)
diff_int = m_bdiff(i,k);
m_Lmatrix(k+n2,i+n2) = 0.0;
sum += x[k]/diff_int;
if (k != i) sum += x[k]*m_astar(i,k)*constant2 /
(m_mw[k]*diff_int);
}
// see Eq. (12.130)
m_Lmatrix(i+n2,i+n2) =
- eightoverpi*m_mw[i]*x[i]*x[i]*m_crot[i] /
(m_cinternal[i]*m_cinternal[i]*GasConstant*m_visc[i]*m_rotrelax[i])
- constant1*sum;
} else {
for (size_t k = 0; k < m_nsp; k++) {
m_Lmatrix(i+n2,i+n2) = 1.0;
}
}
}
}
}
#endif

View file

@ -10,7 +10,6 @@
#include "cantera/transport/MultiTransport.h"
#include "cantera/numerics/ctlapack.h"
#include "cantera/base/utilities.h"
#include "L_matrix.h"
#include "cantera/transport/TransportParams.h"
#include "cantera/thermo/IdealGasPhase.h"
#include "cantera/transport/TransportFactory.h"
@ -583,4 +582,200 @@ void MultiTransport::updateThermal_T()
m_thermal_tlast = m_thermo->temperature();
}
//! Constant to compare dimensionless heat capacities against zero
static const doublereal Min_C_Internal = 0.001;
bool MultiTransport::hasInternalModes(size_t j)
{
return (m_cinternal[j] > Min_C_Internal);
}
void MultiTransport::eval_L0000(const doublereal* const x)
{
doublereal prefactor = 16.0*m_temp/25.0;
doublereal sum;
for (size_t i = 0; i < m_nsp; i++) {
// subtract-off the k=i term to account for the first delta
// function in Eq. (12.121)
sum = -x[i]/m_bdiff(i,i);
for (size_t k = 0; k < m_nsp; k++) {
sum += x[k]/m_bdiff(i,k);
}
sum /= m_mw[i];
for (size_t j = 0; j != m_nsp; ++j) {
m_Lmatrix(i,j) = prefactor * x[j]
* (m_mw[j] * sum + x[i]/m_bdiff(i,j));
}
// diagonal term is zero
m_Lmatrix(i,i) = 0.0;
}
}
void MultiTransport::eval_L0010(const doublereal* const x)
{
doublereal prefactor = 1.6*m_temp;
doublereal sum, wj, xj;
for (size_t j = 0; j < m_nsp; j++) {
//constant = prefactor * x[j];
xj = x[j];
wj = m_mw[j];
sum = 0.0;
for (size_t i = 0; i < m_nsp; i++) {
m_Lmatrix(i,j + m_nsp) = - prefactor * x[i] * xj * m_mw[i] *
(1.2 * m_cstar(j,i) - 1.0) /
((wj + m_mw[i]) * m_bdiff(j,i));
// the next term is independent of "j";
// need to do it for the "j,j" term
sum -= m_Lmatrix(i,j+m_nsp);
}
m_Lmatrix(j,j+m_nsp) += sum;
}
}
void MultiTransport::eval_L1000()
{
for (size_t j = 0; j < m_nsp; j++) {
for (size_t i = 0; i < m_nsp; i++) {
m_Lmatrix(i+m_nsp,j) = m_Lmatrix(j,i+m_nsp);
}
}
}
void MultiTransport::eval_L1010(const doublereal* x)
{
const doublereal fiveover3pi = 5.0/(3.0*Pi);
doublereal prefactor = (16.0*m_temp)/25.0;
doublereal constant1, wjsq, constant2, constant3, constant4,
fourmj, threemjsq, sum, sumwij;;
doublereal term1, term2;
for (size_t j = 0; j < m_nsp; j++) {
// get constant terms that depend on just species "j"
constant1 = prefactor*x[j];
wjsq = m_mw[j]*m_mw[j];
constant2 = 13.75*wjsq;
constant3 = m_crot[j]/m_rotrelax[j];
constant4 = 7.5*wjsq;
fourmj = 4.0*m_mw[j];
threemjsq = 3.0*m_mw[j]*m_mw[j];
sum = 0.0;
for (size_t i = 0; i < m_nsp; i++) {
sumwij = m_mw[i] + m_mw[j];
term1 = m_bdiff(i,j) * sumwij*sumwij;
term2 = fourmj*m_astar(i,j)*(1.0 + fiveover3pi*
(constant3 +
(m_crot[i]/m_rotrelax[i]))); // see Eq. (12.125)
m_Lmatrix(i+m_nsp,j+m_nsp) = constant1*x[i]*m_mw[i] /(m_mw[j]*term1) *
(constant2 - threemjsq*m_bstar(i,j)
- term2*m_mw[j]);
sum += x[i] /(term1) *
(constant4 + m_mw[i]*m_mw[i]*
(6.25 - 3.0*m_bstar(i,j)) + term2*m_mw[i]);
}
m_Lmatrix(j+m_nsp,j+m_nsp) -= sum*constant1;
}
}
void MultiTransport::eval_L1001(const doublereal* x)
{
doublereal prefactor = 32.00*m_temp/(5.00*Pi);
doublereal constant, sum;
size_t n2 = 2*m_nsp;
int npoly = 0;
for (size_t j = 0; j < m_nsp; j++) {
// collect terms that depend only on "j"
if (hasInternalModes(j)) {
constant = prefactor*m_mw[j]*x[j]*m_crot[j]/(m_cinternal[j]*m_rotrelax[j]);
sum = 0.0;
for (size_t i = 0; i < m_nsp; i++) {
// see Eq. (12.127)
m_Lmatrix(i+m_nsp,j+n2) = constant * m_astar(j,i) * x[i] /
((m_mw[j] + m_mw[i]) * m_bdiff(j,i));
sum += m_Lmatrix(i+m_nsp,j+n2);
}
npoly++;
m_Lmatrix(j+m_nsp,j+n2) += sum;
} else {
for (size_t i = 0; i < m_nsp; i++) {
m_Lmatrix(i+m_nsp,j+n2) = 0.0;
}
}
}
}
void MultiTransport::eval_L0001()
{
size_t n2 = 2*m_nsp;
for (size_t j = 0; j < m_nsp; j++) {
for (size_t i = 0; i < m_nsp; i++) {
m_Lmatrix(i,j+n2) = 0.0;
}
}
}
void MultiTransport::eval_L0100()
{
size_t n2 = 2*m_nsp;
for (size_t j = 0; j < m_nsp; j++)
for (size_t i = 0; i < m_nsp; i++) {
m_Lmatrix(i+n2,j) = 0.0; // see Eq. (12.123)
}
}
void MultiTransport::eval_L0110()
{
size_t n2 = 2*m_nsp;
for (size_t j = 0; j < m_nsp; j++)
for (size_t i = 0; i < m_nsp; i++) {
m_Lmatrix(i+n2,j+m_nsp) = m_Lmatrix(j+m_nsp,i+n2); // see Eq. (12.123)
}
}
void MultiTransport::eval_L0101(const doublereal* x)
{
const doublereal fivepi = 5.00*Pi;
const doublereal eightoverpi = 8.0 / Pi;
doublereal prefactor = 4.00*m_temp;
size_t n2 = 2*m_nsp;
doublereal constant1, constant2, diff_int, sum;
for (size_t i = 0; i < m_nsp; i++) {
if (hasInternalModes(i)) {
// collect terms that depend only on "i"
constant1 = prefactor*x[i]/m_cinternal[i];
constant2 = 12.00*m_mw[i]*m_crot[i] /
(fivepi*m_cinternal[i]*m_rotrelax[i]);
sum = 0.0;
for (size_t k = 0; k < m_nsp; k++) {
// see Eq. (12.131)
diff_int = m_bdiff(i,k);
m_Lmatrix(k+n2,i+n2) = 0.0;
sum += x[k]/diff_int;
if (k != i) sum += x[k]*m_astar(i,k)*constant2 /
(m_mw[k]*diff_int);
}
// see Eq. (12.130)
m_Lmatrix(i+n2,i+n2) =
- eightoverpi*m_mw[i]*x[i]*x[i]*m_crot[i] /
(m_cinternal[i]*m_cinternal[i]*GasConstant*m_visc[i]*m_rotrelax[i])
- constant1*sum;
} else {
for (size_t k = 0; k < m_nsp; k++) {
m_Lmatrix(i+n2,i+n2) = 1.0;
}
}
}
}
}