36 namespace laminarFlameSpeedModels
52 Foam::laminarFlameSpeedModels::SCOPE::polynomial::polynomial
54 const dictionary& polyDict
57 FixedList<scalar, 7>(polyDict.lookup(
"coefficients")),
58 ll(polyDict.
get<scalar>(
"lowerLimit")),
59 ul(polyDict.
get<scalar>(
"upperLimit")),
60 llv(polyPhi(ll, *this)),
61 ulv(polyPhi(ul, *this)),
66 Foam::laminarFlameSpeedModels::SCOPE::SCOPE
68 const dictionary&
dict,
69 const psiuReactionThermo& ct
72 laminarFlameSpeed(
dict, ct),
82 ).optionalSubDict(typeName +
"Coeffs")
86 coeffsDict_.getCompat<scalar>
88 "lowerFlammabilityLimit",
89 {{
"lowerFlamabilityLimit", 1712}}
94 coeffsDict_.getCompat<scalar>
96 "upperFlammabilityLimit",
97 {{
"upperFlamabilityLimit", 1712}}
100 SuPolyL_(coeffsDict_.subDict(
"lowerSuPolynomial")),
101 SuPolyU_(coeffsDict_.subDict(
"upperSuPolynomial")),
102 Texp_(coeffsDict_.get<scalar>(
"Texp")),
103 pexp_(coeffsDict_.get<scalar>(
"pexp")),
104 MaPolyL_(coeffsDict_.subDict(
"lowerMaPolynomial")),
105 MaPolyU_(coeffsDict_.subDict(
"upperMaPolynomial"))
107 SuPolyL_.ll =
max(SuPolyL_.ll, LFL_) + SMALL;
108 SuPolyU_.ul =
min(SuPolyU_.ul, UFL_) - SMALL;
110 SuPolyL_.lu = 0.5*(SuPolyL_.ul + SuPolyU_.ll);
111 SuPolyU_.lu = SuPolyL_.lu - SMALL;
113 MaPolyL_.lu = 0.5*(MaPolyL_.ul + MaPolyU_.ll);
114 MaPolyU_.lu = MaPolyL_.lu - SMALL;
118 Info<<
"phi Su (T = Tref, p = pref)" <<
endl;
120 for (
int i=0; i<
n; i++)
122 scalar
phi = (2.0*i)/
n;
137 inline Foam::scalar Foam::laminarFlameSpeedModels::SCOPE::polyPhi
143 scalar
x =
phi - 1.0;
149 +
x*(a[1] +
x*(a[2] +
x*(a[3] +
x*(a[4] +
x*(a[5] +
x*a[6])))))
154 inline Foam::scalar Foam::laminarFlameSpeedModels::SCOPE::SuRef
159 if (phi < LFL_ || phi > UFL_)
164 else if (
phi < SuPolyL_.ll)
168 return SuPolyL_.llv*(
phi - LFL_)/(SuPolyL_.ll - LFL_);
170 else if (
phi > SuPolyU_.ul)
174 return SuPolyU_.ulv*(UFL_ -
phi)/(UFL_ - SuPolyU_.ul);
176 else if (
phi < SuPolyL_.lu)
179 return polyPhi(
phi, SuPolyL_);
181 else if (
phi > SuPolyU_.lu)
184 return polyPhi(
phi, SuPolyU_);
190 <<
" cannot be handled by SCOPE function with the " 204 if (
phi < MaPolyL_.ll)
209 else if (
phi > MaPolyU_.ul)
214 else if (
phi < SuPolyL_.lu)
217 return polyPhi(
phi, MaPolyL_);
219 else if (
phi > SuPolyU_.lu)
222 return polyPhi(
phi, MaPolyU_);
228 <<
" cannot be handled by SCOPE function with the " 237 inline Foam::scalar Foam::laminarFlameSpeedModels::SCOPE::Su0pTphi
244 static const scalar Tref = 300.0;
245 static const scalar pRef = 1.013e5;
247 return SuRef(
phi)*
pow((Tu/Tref), Texp_)*
pow((
p/pRef), pexp_);
258 tmp<volScalarField> tSu0
279 Su0[celli] = Su0pTphi(
p[celli], Tu[celli],
phi);
288 const scalarField& Tup = Tu.boundaryField()[patchi];
292 Su0p[facei] = Su0pTphi(
pp[facei], Tup[facei],
phi);
307 tmp<volScalarField> tSu0
328 Su0[celli] = Su0pTphi(
p[celli], Tu[celli],
phi[celli]);
337 const scalarField& Tup = Tu.boundaryField()[patchi];
361 tmp<volScalarField> tMa
368 phi.time().timeName(),
382 ma[celli] = Ma(
phi[celli]);
394 map[facei] = Ma(phip[facei]);
405 if (psiuReactionThermo_.composition().contains(
"ft"))
407 const volScalarField& ft = psiuReactionThermo_.composition().Y(
"ft");
413 "stoichiometricAirFuelMassRatio",
dimless, psiuReactionThermo_
414 )*ft/(scalar(1) - ft)
419 const fvMesh&
mesh = psiuReactionThermo_.p().mesh();
421 return tmp<volScalarField>
428 mesh.time().timeName(),
444 if (psiuReactionThermo_.composition().contains(
"ft"))
446 const volScalarField& ft = psiuReactionThermo_.composition().Y(
"ft");
450 psiuReactionThermo_.p(),
451 psiuReactionThermo_.Tu(),
454 "stoichiometricAirFuelMassRatio",
dimless, psiuReactionThermo_
455 )*ft/(scalar(1) - ft)
462 psiuReactionThermo_.p(),
463 psiuReactionThermo_.Tu(),
List< ReturnType > get(const UPtrList< T > &list, const AccessOp &aop)
List of values generated by applying the access operation to each list item.
errorManipArg< error, int > exit(error &err, const int errNo=1)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Ostream & endl(Ostream &os)
Add newline and flush stream.
tmp< volScalarField > operator()() const
Return the laminar flame speed [m/s].
GeometricBoundaryField< scalar, fvPatchField, volMesh > Boundary
Type of boundary fields.
Ignore writing from objectRegistry::writeObject()
const dimensionSet dimless
Dimensionless.
Macros for easy insertion into run-time selection tables.
#define forAll(list, i)
Loop across all elements in list.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
addToRunTimeSelectionTable(laminarFlameSpeed, constant, dictionary)
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
int debug
Static debugging option.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
messageStream Info
Information stream (stdout output on master, null elsewhere)
Internal & ref(const bool updateAccessTime=true)
Same as internalFieldRef()
A class for managing temporary objects.
tmp< volScalarField > Ma() const
Return the Markstein number.
defineTypeNameAndDebug(constant, 0)
uindirectPrimitivePatch pp(UIndirectList< face >(mesh.faces(), faceLabels), mesh.points())
static constexpr const zero Zero
Global zero (0)
const dimensionSet dimVelocity