34 namespace laminarFlameSpeedModels
50 Foam::laminarFlameSpeedModels::SCOPE::polynomial::polynomial
52 const dictionary& polyDict
55 FixedList<scalar, 7>(polyDict.lookup(
"coefficients")),
56 ll(polyDict.lookup<scalar>(
"lowerLimit")),
57 ul(polyDict.lookup<scalar>(
"upperLimit")),
58 llv(polyPhi(ll, *this)),
59 ulv(polyPhi(ul, *this)),
80 dict.lookup(
"fuelFile")
83 ).optionalSubDict(typeName +
"Coeffs")
85 LFL_(coeffsDict_.lookup<scalar>(
"lowerFlamabilityLimit")),
86 UFL_(coeffsDict_.lookup<scalar>(
"upperFlamabilityLimit")),
87 SuPolyL_(coeffsDict_.subDict(
"lowerSuPolynomial")),
88 SuPolyU_(coeffsDict_.subDict(
"upperSuPolynomial")),
89 Texp_(coeffsDict_.lookup<scalar>(
"Texp")),
90 pexp_(coeffsDict_.lookup<scalar>(
"pexp")),
91 MaPolyL_(coeffsDict_.subDict(
"lowerMaPolynomial")),
92 MaPolyU_(coeffsDict_.subDict(
"upperMaPolynomial"))
94 SuPolyL_.ll =
max(SuPolyL_.ll, LFL_) + small;
95 SuPolyU_.ul =
min(SuPolyU_.ul, UFL_) - small;
97 SuPolyL_.lu = 0.5*(SuPolyL_.ul + SuPolyU_.ll);
98 SuPolyU_.lu = SuPolyL_.lu - small;
100 MaPolyL_.lu = 0.5*(MaPolyL_.ul + MaPolyU_.ll);
101 MaPolyU_.lu = MaPolyL_.lu - small;
105 Info<<
"phi Su (T = Tref, p = pref)" <<
endl;
107 for (
int i=0; i<
n; i++)
109 scalar phi = (2.0*i)/
n;
124 inline Foam::scalar Foam::laminarFlameSpeedModels::SCOPE::polyPhi
130 scalar
x = phi - 1.0;
136 +
x*(a[1] +
x*(a[2] +
x*(a[3] +
x*(a[4] +
x*(a[5] +
x*a[6])))))
141 inline Foam::scalar Foam::laminarFlameSpeedModels::SCOPE::SuRef
146 if (phi < LFL_ || phi > UFL_)
151 else if (phi < SuPolyL_.ll)
155 return SuPolyL_.llv*(phi - LFL_)/(SuPolyL_.ll - LFL_);
157 else if (phi > SuPolyU_.ul)
161 return SuPolyU_.ulv*(UFL_ - phi)/(UFL_ - SuPolyU_.ul);
163 else if (phi < SuPolyL_.lu)
166 return polyPhi(phi, SuPolyL_);
168 else if (phi > SuPolyU_.lu)
171 return polyPhi(phi, SuPolyU_);
177 <<
" cannot be handled by SCOPE function with the "
191 if (phi < MaPolyL_.ll)
196 else if (phi > MaPolyU_.ul)
201 else if (phi < SuPolyL_.lu)
204 return polyPhi(phi, MaPolyL_);
206 else if (phi > SuPolyU_.lu)
209 return polyPhi(phi, MaPolyU_);
215 <<
" cannot be handled by SCOPE function with the "
224 inline Foam::scalar Foam::laminarFlameSpeedModels::SCOPE::Su0pTphi
231 static const scalar Tref = 300.0;
232 static const scalar pRef = 1.013e5;
234 return SuRef(phi)*
pow((Tu/Tref), Texp_)*
pow((
p/pRef), pexp_);
245 tmp<volScalarField> tSu0
259 Su0[celli] = Su0pTphi(
p[celli], Tu[celli], phi);
272 Su0p[facei] = Su0pTphi(pp[facei], Tup[facei], phi);
287 tmp<volScalarField> tSu0
301 Su0[celli] = Su0pTphi(
p[celli], Tu[celli], phi[celli]);
334 tmp<volScalarField> tMa
348 ma[celli] = Ma(phi[celli]);
360 map[facei] = Ma(phip[facei]);
371 if (psiuMulticomponentThermo_.containsSpecie(
"ft"))
379 "stoichiometricAirFuelMassRatio",
381 psiuMulticomponentThermo_.properties()
382 )*ft/(scalar(1) - ft)
387 const fvMesh& mesh = psiuMulticomponentThermo_.p().
mesh();
405 if (psiuMulticomponentThermo_.containsSpecie(
"ft"))
411 psiuMulticomponentThermo_.p(),
412 psiuMulticomponentThermo_.Tu(),
415 "stoichiometricAirFuelMassRatio",
417 psiuMulticomponentThermo_.properties()
418 )*ft/(scalar(1) - ft)
425 psiuMulticomponentThermo_.p(),
426 psiuMulticomponentThermo_.Tu(),
#define forAll(list, i)
Loop across all elements in list.
Macros for easy insertion into run-time selection tables.
Generic GeometricField class.
GeometricBoundaryField< Type, PatchField, GeoMesh > Boundary
Type of the boundary field.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
A list of keyword definitions, which are a keyword followed by any number of values (e....
A class for handling file names.
Mesh data needed to do the Finite Volume discretisation.
const polyMesh & mesh() const
Return reference to polyMesh.
Laminar flame speed obtained from the SCOPE correlation.
tmp< volScalarField > Ma() const
Return the Markstein number.
tmp< volScalarField > operator()() const
Return the laminar flame speed [m/s].
SCOPE(const dictionary &, const psiuMulticomponentThermo &)
Construct from dictionary and psiuMulticomponentThermo.
Abstract class for laminar flame speed.
Base-class for combustion fluid thermodynamic properties based on compressibility.
A class for managing temporary objects.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
addToRunTimeSelectionTable(laminarFlameSpeed, constant, dictionary)
defineTypeNameAndDebug(constant, 0)
errorManipArg< error, int > exit(error &err, const int errNo=1)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Ostream & endl(Ostream &os)
Add newline and flush stream.
const dimensionSet dimless
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
VolField< scalar > volScalarField
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
const dimensionSet dimVelocity
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.