interpolate SLFM library to integral grid and removed redundant pdf calculation
This commit is contained in:
ignis
parent
44698c3f7a
commit
952de05fd0
3 changed files with 89 additions and 27 deletions
|
|
@ -3,15 +3,27 @@ class BetaPDF
|
|||
|
||||
//private variable
|
||||
|
||||
//raw parameter mixture fraction and mf variance
|
||||
scalar mf_, mfVar_;
|
||||
|
||||
//beta-pdf parameter alpha and beta
|
||||
scalar alpha_, beta_;
|
||||
|
||||
//cutting point and number of eta-space
|
||||
scalarField etaCut_, N_;
|
||||
//cutting points
|
||||
scalarField etaCut_;
|
||||
|
||||
//number of eta-space
|
||||
labelList N_;
|
||||
|
||||
//detailed integration space
|
||||
scalarField etaSpace_;
|
||||
|
||||
//pdf numerators
|
||||
scalarField pdfNum_;
|
||||
|
||||
//pdf denominator
|
||||
scalar pdfDen_;
|
||||
|
||||
//detailed integration space, part
|
||||
typedef List<scalarField> scalarFieldArray1d;
|
||||
scalarFieldArray1d etaPart_;
|
||||
|
|
@ -59,6 +71,9 @@ public:
|
|||
}
|
||||
|
||||
AMCfine_ = AMC(etaSpace_);
|
||||
|
||||
pdfNum_ = scalarField(etaSpace_.size(), 0.0);
|
||||
pdfDen_ = 1.0;
|
||||
}
|
||||
|
||||
// Destructor
|
||||
|
|
@ -67,8 +82,11 @@ public:
|
|||
// Member functions
|
||||
|
||||
// set beta-pdf parameter alpha_ and beta_
|
||||
void setParameter(scalar& mf, scalar& mfVar)
|
||||
void setParameter(const scalar mf, const scalar mfVar)
|
||||
{
|
||||
mf_ = mf;
|
||||
mfVar_ = mfVar;
|
||||
|
||||
fdelta_ = false;
|
||||
delta_ox = false;
|
||||
delta_fu = false;
|
||||
|
|
@ -92,10 +110,17 @@ public:
|
|||
{
|
||||
delta_fu = true;
|
||||
}
|
||||
|
||||
limit_ab();
|
||||
}
|
||||
|
||||
pdfNum_ = etaFunc(alpha_, beta_, etaSpace_);
|
||||
pdfDen_ = integrate(pdfNum_)
|
||||
+ Foam::pow(etaSpace_[0], alpha_)/(alpha_ + SMALL)
|
||||
+ Foam::pow(etaSpace_[0], beta_)/(beta_ + SMALL);
|
||||
}
|
||||
// beta-pdf weighted integration for given mf, mfVar and f
|
||||
scalar evaluate(scalar& mf, scalar& mfVar, scalarField& etaValue, scalarField& f)
|
||||
scalar evaluate(scalar mf, scalar mfVar, scalarField& etaValue, scalarField& f)
|
||||
{
|
||||
scalar result(0);
|
||||
|
||||
|
|
@ -103,8 +128,6 @@ public:
|
|||
|
||||
if(fdelta_ == false)
|
||||
{
|
||||
limit_ab();
|
||||
|
||||
scalar num(0), den(0);
|
||||
for(label i=0 ; i<etaCut_.size()-1 ; i++)
|
||||
{
|
||||
|
|
@ -130,7 +153,45 @@ public:
|
|||
|
||||
return result;
|
||||
}
|
||||
scalar value(scalar& mf, scalar& mfVar, scalar& etaValue)
|
||||
|
||||
// beta-pdf weighted integration for given mf, mfVar and f
|
||||
scalar betaIntegrate(const scalarField& f)
|
||||
{
|
||||
scalar result(0);
|
||||
|
||||
if(fdelta_)
|
||||
{
|
||||
result = interpolateXY(mf_, etaSpace_, f);
|
||||
}
|
||||
else
|
||||
{
|
||||
scalar num = integrate (f * pdfNum_)
|
||||
+ f.first()*Foam::pow(etaSpace_[0], alpha_)/(alpha_ + SMALL)
|
||||
+ f.last()*Foam::pow(etaSpace_[0], beta_)/(beta_ + SMALL);
|
||||
result = num/pdfDen_;
|
||||
}
|
||||
|
||||
return result;
|
||||
}
|
||||
|
||||
scalar integrate(const scalarField& f)
|
||||
{
|
||||
scalar total = 0.0;
|
||||
|
||||
forAll(N_, i)
|
||||
{
|
||||
const labelList::subList prev(N_, i);
|
||||
const label baseI = sum(prev);
|
||||
|
||||
const scalarField::subField subInterval(f, N_[i]+1, baseI);
|
||||
|
||||
total += simps(etaCut_[i], etaCut_[i+1], N_[i], subInterval);
|
||||
}
|
||||
|
||||
return total;
|
||||
}
|
||||
|
||||
scalar value(scalar mf, scalar mfVar, scalar etaValue)
|
||||
{
|
||||
scalar result(0);
|
||||
|
||||
|
|
@ -138,8 +199,6 @@ public:
|
|||
|
||||
if(fdelta_ == false)
|
||||
{
|
||||
limit_ab();
|
||||
|
||||
scalar den(0);
|
||||
for(label i=0 ; i<etaCut_.size()-1 ; i++)
|
||||
{
|
||||
|
|
@ -172,21 +231,21 @@ public:
|
|||
alpha_ = (1.0+fmax*(beta_-2.0))/(1.0-fmax);
|
||||
}
|
||||
}
|
||||
scalarField etaFunc(scalar& a, scalar& b, scalarField& eta)
|
||||
scalarField etaFunc(scalar a, scalar b, scalarField& eta)
|
||||
{
|
||||
return pow(eta, a-1.0)*pow(1.0-eta, b-1.0);
|
||||
}
|
||||
scalar etaFunc(scalar& a, scalar& b, scalar& eta)
|
||||
scalar etaFunc(scalar a, scalar b, scalar eta)
|
||||
{
|
||||
return Foam::pow(eta, a-1.0)*Foam::pow(1.0-eta, b-1.0);
|
||||
}
|
||||
//extended Simpson's rule (Numerical recipes, 2nd Ed. p.128)
|
||||
//for equally spaced and even N intervals (or odd N+1 points)
|
||||
scalar simps(scalar& xl, scalar& xh, scalar& N, scalarField& fx)
|
||||
scalar simps(const scalar xl, const scalar xh, const label N, const UList<scalar>& fx)
|
||||
{
|
||||
scalar evensum(0.0), oddsum(0.0), sum(0.0);
|
||||
|
||||
scalar h = (xh - xl)/N;
|
||||
scalar h = (xh - xl)/scalar(N);
|
||||
|
||||
for(label i=0 ; i<fx.size() ; i++)
|
||||
{
|
||||
|
|
@ -228,11 +287,11 @@ public:
|
|||
|
||||
return result;
|
||||
}
|
||||
scalar AMC(scalar& eta)
|
||||
scalar AMC(scalar eta)
|
||||
{
|
||||
return interpolateXY(eta, etaSpace_, AMCfine_);
|
||||
}
|
||||
void C1coeff(scalar& mf, scalarField& varValue, scalarField& C1table)
|
||||
void C1coeff(scalar mf, scalarField& varValue, scalarField& C1table)
|
||||
{
|
||||
scalar maxVar = mf*(1.0-mf);
|
||||
scalarField x, fx;
|
||||
|
|
|
|||
|
|
@ -47,18 +47,21 @@ for(label j=1 ; j<etamax-1; j++)
|
|||
Info<<"Construct Temp and Rgas table"<<endl;
|
||||
for(label j=0 ; j<etamax ; j++)
|
||||
{
|
||||
scalar mf = etaValue[j];
|
||||
const scalar mf = etaValue[j];
|
||||
for(label v=0 ; v<=NVar ; v++)
|
||||
{
|
||||
scalar mfVar = varValue[v]*mf*(1.0-mf);
|
||||
const scalar mfVar = varValue[v]*mf*(1.0-mf);
|
||||
|
||||
bpdf.setParameter(mf, mfVar);
|
||||
|
||||
for(label n=0 ; n<NstList.size() ; n++)
|
||||
{
|
||||
Ttable[j][v][n] = bpdf.evaluate(mf, mfVar, etaValue, T_SLFM[n]);
|
||||
Rtable[j][v][n] = bpdf.evaluate(mf, mfVar, etaValue, Rgas_SLFM[n]);
|
||||
Ttable[j][v][n] = bpdf.betaIntegrate(T_SLFM[n]);
|
||||
Rtable[j][v][n] = bpdf.betaIntegrate(Rgas_SLFM[n]);
|
||||
forAll(postSpecieIndices, yi)
|
||||
{
|
||||
const label y = postSpecieIndices[yi];
|
||||
Ytable[yi][j][v][n] = bpdf.evaluate(mf, mfVar, etaValue, Y_SLFM[y][n]);
|
||||
Ytable[yi][j][v][n] = bpdf.betaIntegrate(Y_SLFM[y][n]);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
|
|
|||
|
|
@ -12,7 +12,7 @@ typedef List<scalarField> scalarFieldArray1d;
|
|||
typedef List<scalarFieldArray1d> scalarFieldArray2d;
|
||||
typedef List<scalarFieldArray2d> scalarFieldArray3d;
|
||||
|
||||
CMC::FlameStructure slfmLibrary(SLFMdict, etaValue.size(), Y.size());
|
||||
CMC::FlameStructure slfmLibrary(SLFMdict, etaGrid.size(), Y.size());
|
||||
|
||||
scalarFieldArray2d& Y_SLFM (slfmLibrary.Y());
|
||||
scalarFieldArray2d& W_SLFM (slfmLibrary.W());
|
||||
|
|
@ -88,7 +88,7 @@ for(label n=0 ; n<NstList.size() ; n++)
|
|||
Y_temp[j] = max(0.0, Yj);
|
||||
}
|
||||
|
||||
Y_SLFM[spI][n] = interpolateXY(etaValue, etaValue_SLFM, Y_temp);
|
||||
Y_SLFM[spI][n] = interpolateXY(etaGrid, etaValue_SLFM, Y_temp);
|
||||
|
||||
//read reaction rate (from SLFM library)
|
||||
for(label j=0 ; j<etamax_SLFM ; j++)
|
||||
|
|
@ -96,7 +96,7 @@ for(label n=0 ; n<NstList.size() ; n++)
|
|||
fin>>W_temp[j];
|
||||
}
|
||||
|
||||
W_SLFM[spI][n] = interpolateXY(etaValue, etaValue_SLFM, W_temp);
|
||||
W_SLFM[spI][n] = interpolateXY(etaGrid, etaValue_SLFM, W_temp);
|
||||
|
||||
//blank line
|
||||
fin.getLine(gbg);
|
||||
|
|
@ -115,7 +115,7 @@ for(label n=0 ; n<NstList.size() ; n++)
|
|||
Y_temp[j] = max(0.0, Yj);
|
||||
}
|
||||
|
||||
T_SLFM[n] = interpolateXY(etaValue, etaValue_SLFM, Y_temp);
|
||||
T_SLFM[n] = interpolateXY(etaGrid, etaValue_SLFM, Y_temp);
|
||||
|
||||
//read heat source
|
||||
for(label j=0 ; j<etamax_SLFM ; j++)
|
||||
|
|
@ -123,7 +123,7 @@ for(label n=0 ; n<NstList.size() ; n++)
|
|||
fin>>W_temp[j];
|
||||
}
|
||||
|
||||
Q_SLFM[n] = interpolateXY(etaValue, etaValue_SLFM, W_temp);
|
||||
Q_SLFM[n] = interpolateXY(etaGrid, etaValue_SLFM, W_temp);
|
||||
|
||||
//blank line
|
||||
fin.getLine(gbg);
|
||||
|
|
@ -140,7 +140,7 @@ for(label n=0 ; n<NstList.size() ; n++)
|
|||
}
|
||||
//Info<<rho_temp<<endl;
|
||||
|
||||
rho_SLFM[n] = interpolateXY(etaValue, etaValue_SLFM, rho_temp);
|
||||
rho_SLFM[n] = interpolateXY(etaGrid, etaValue_SLFM, rho_temp);
|
||||
|
||||
//blank line
|
||||
fin.getLine(gbg);
|
||||
|
|
@ -155,7 +155,7 @@ for(label n=0 ; n<NstList.size() ; n++)
|
|||
fin>>h_temp[j];
|
||||
}
|
||||
|
||||
h_SLFM[n] = interpolateXY(etaValue, etaValue_SLFM, h_temp);
|
||||
h_SLFM[n] = interpolateXY(etaGrid, etaValue_SLFM, h_temp);
|
||||
|
||||
for(label i=0 ; i<Y.size() ; i++)
|
||||
{
|
||||
|
|
|
|||
Loading…
Add table
Reference in a new issue