replaced hard coded constants

This commit is contained in:
ignis 2018-02-01 23:10:09 +09:00
parent 09c95558db
commit f752bb2a50
3 changed files with 36 additions and 51 deletions

View file

@ -26,6 +26,7 @@ License
#include "MossBrookesSoot.H"
#include "addToRunTimeSelectionTable.H"
#include "fvm.H"
#include "constants.H"
#include "basicSpecieMixture.H"
@ -269,6 +270,7 @@ Foam::radiation::MossBrookesSoot::MossBrookesSoot
mesh_,
dimensionedScalar("SootMeanDiameter", dimLength, 0.0)
),
Nnorm("Nnorm", dimless, 1e15),
PrtSoot(readScalar(coeffsDict_.lookup("TurbulentPrandtlNumberSoot"))),
rhoSoot(readScalar(coeffsDict_.lookup("SootDensity"))),
Calpha(readScalar(coeffsDict_.lookup("Calpha"))),
@ -279,9 +281,7 @@ Foam::radiation::MossBrookesSoot::MossBrookesSoot
etaColl(readScalar(coeffsDict_.lookup("etaColl"))),
Talpha(readScalar(coeffsDict_.lookup("Talpha"))),
Tgamma(readScalar(coeffsDict_.lookup("Tgamma"))),
Navog(readScalar(coeffsDict_.lookup("AvogadroNumber"))),
Cbeta(readScalar(coeffsDict_.lookup("Cbeta"))),
Boltzmann(readScalar(coeffsDict_.lookup("BoltzmannConst"))),
radiationActive_(readBool(coeffsDict_.lookup("Radiationeffect"))),
cliqueNuc(readScalar(coeffsDict_.lookup("NucExpCorrectConstant"))),
cliqueMassGrow(readScalar(coeffsDict_.lookup("MassGrowCorrectConstant")))
@ -381,34 +381,13 @@ void Foam::radiation::MossBrookesSoot::correct()
NucConcEqn.solve();
scalar minusindicator1 = 0.0;
scalar minusindicator2 = 0.0;
forAll(rho_, celli)
{
if (soot_[celli] < 0.0) {
soot_[celli] = 0.0;
minusindicator1 = 1.0;
}
if (NucConc_[celli] < 0.0) {
NucConc_[celli] =0.0;
minusindicator2 = 1.0;
}
soot_.max(0.0);
NucConc_.max(0.0);
NumDen_[celli] = rho_[celli]*(NucConc_[celli]*1e+15);
NumDen_ = rho_*NucConc_*Nnorm;
sootVF_[celli] = soot_[celli]*rho_[celli]/rhoSoot;
}
if (minusindicator1 == 1.0 and minusindicator2 == 0.0) {
Info<< "Minus soot value! corrected to zero" <<endl;
}
else if (minusindicator1 == 0.0 and minusindicator2 ==1.0) {
Info<< "Minus NucConc value! corrected to zero" <<endl;
}
else if (minusindicator1 == 1.0 and minusindicator2 == 1.0) {
Info<< "Minus soot and NucConc value! corrected to zero" << endl;
}
sootVF_ = soot_*rho_/dimensionedScalar("rhoSoot", dimDensity, rhoSoot);
}
@ -430,6 +409,12 @@ void Foam::radiation::MossBrookesSoot::calcSource()
const volScalarField XSgs_ (mixture_.W()*Ygrow_/mixture_.W(growI));
const volScalarField XOxid_(mixture_.W()*Yoxid_/mixture_.W(oxidI));
const dimensionedScalar NA (Foam::constant::physicoChemical::NA);
const dimensionedScalar kB (Foam::constant::physicoChemical::k);
const dimensionedScalar R (Foam::constant::physicoChemical::R);
const scalar pi314 = Foam::constant::mathematical::pi;
forAll(rho_, celli)
{
const scalar XPreci = XPrec_[celli];
@ -441,36 +426,38 @@ void Foam::radiation::MossBrookesSoot::calcSource()
const scalar sooti = soot_[celli];
const scalar NumDeni = NumDen_[celli];
const scalar MolarConcPreci = (XPreci*pi)/(8314.5*Ti);
const scalar MolarConcSgsi = (XSgsi*pi)/(8314.5*Ti);
const scalar MolarConcOxidi = (XOxidi*pi)/(8314.5*Ti);
// concentrations in [kmol/m3]
const scalar MolarConcPreci = (XPreci*pi)/(1e3*R.value()*Ti);
const scalar MolarConcSgsi = (XSgsi *pi)/(1e3*R.value()*Ti);
const scalar MolarConcOxidi = (XOxidi*pi)/(1e3*R.value()*Ti);
scalar SootMeanDiameteri = 0.0;
if (NumDeni == 0.0) {
if (NumDeni == 0.0)
{
sootMeanDiameter_[celli] = 0.0;
Sf_[celli] = 0.0;
}
else {
const scalar parameter = (6.0/3.1416)*(rhoi/rhoSoot)*(sooti/NumDeni);
sootMeanDiameter_[celli] = pow(parameter,0.33333);
SootMeanDiameteri = sootMeanDiameter_[celli];
Sf_[celli] = (3.1416)*pow(SootMeanDiameteri,2.0)*NumDeni;
else
{
const scalar parameter = (6.0/pi314)*(rhoi/rhoSoot)*(sooti/NumDeni);
sootMeanDiameter_[celli] = pow(parameter, 1./3.);
Sf_[celli] = pi314*pow(sootMeanDiameter_[celli],2.0)*NumDeni;
}
const scalar sootMeanDiameteri = sootMeanDiameter_[celli];
Snuc_[celli] = Mp*Calpha*exp(-cliqueNuc*Talpha/Ti)*MolarConcPreci;
Smassgrow_[celli] = Cgamma*MolarConcSgsi*exp(-cliqueMassGrow*Tgamma/Ti)*Sf_[celli];
Soxid_[celli] = Coxid*Cw*etaColl*MolarConcOxidi*pow(Ti,0.5)*Sf_[celli];
Soxid_[celli] = Coxid*Cw*etaColl*MolarConcOxidi*sqrt(Ti)*Sf_[celli];
Snet_[celli] = Snuc_[celli] + Smassgrow_[celli] - Soxid_[celli];
NSnuc_[celli] = (Calpha*Navog*MolarConcPreci*exp(-cliqueNuc*Talpha/Ti))*(1e-15);
NSnuc_[celli] = (Calpha*NA.value()*MolarConcPreci*exp(-cliqueNuc*Talpha/Ti))/Nnorm.value();
aggParam_[celli] = 24.0*Boltzmann*Ti/rhoSoot;
aggParam_[celli] = 24.0*kB.value()*Ti/rhoSoot;
NSagg_[celli] = (Cbeta*pow((aggParam_[celli]),0.5)*pow(SootMeanDiameteri,0.5)*pow(NumDeni,2.0))*(1e-15);
NSagg_[celli] = (Cbeta*sqrt(aggParam_[celli])*sqrt(sootMeanDiameteri)*pow(NumDeni,2.0))/Nnorm.value();
NSnet_[celli] = NSnuc_[celli] - NSagg_[celli];
}

View file

@ -159,6 +159,9 @@ class MossBrookesSoot
//- Density of soot particle
scalar rhoSoot;
//- soot number density normalization factor
const dimensionedScalar Nnorm;
//- Model constant for soot inception rate
scalar Calpha;
@ -183,15 +186,9 @@ class MossBrookesSoot
//- Activation temperature for soot growth rate
scalar Tgamma;
//- Avogadro's number
scalar Navog;
//- Model constant for coagulation rate
scalar Cbeta;
//- Boltzmann constant
scalar Boltzmann;
//- Soot radiative effect on/off flag
bool radiationActive_;

View file

@ -348,12 +348,13 @@ void Foam::radiation::khanGreeveSoot::calcSource()
S1_[celli] = A_*rhoi*sooti*epsiloni/ki;
if (sooti == 0.0 and Yfueli == 0.0 ){
if (sooti == 0.0 and Yfueli == 0.0 )
{
S2_[celli] = 0.0;
}
else {
else
{
S2_[celli] = A_*rhoi*(Yoxi/nuSoot_)*(sooti*nuSoot_/(sooti*nuSoot_+(Yfueli)*nuFuel_))*epsiloni/ki;
}
Scomb_[celli] = min(S1_[celli], S2_[celli]);