44 namespace functionObjects
50 BilgerMixtureFraction,
59 void Foam::functionObjects::BilgerMixtureFraction::calcBilgerMixtureFraction()
79 mesh_.objectRegistry::store(resultPtr);
81 auto& f_Bilger = *resultPtr;
83 auto&
Y = thermo_.
Y();
85 f_Bilger = -o2RequiredOx_;
90 *(nAtomsC_[i] + nAtomsS_[i] + 0.25*nAtomsH_[i] - 0.5*nAtomsO_[i])
93 f_Bilger /= o2RequiredFuelOx_;
94 f_Bilger.clamp_range(zero_one{});
98 bool Foam::functionObjects::BilgerMixtureFraction::readComposition
100 const dictionary& subDict,
104 auto&
Y = thermo_.Y();
111 subDict.getCheckOrDefault<scalar>
119 if (
sum(comp) < SMALL)
122 <<
"No composition is given" <<
nl 123 <<
"Valid species are:" <<
nl 130 const word fractionBasisType
132 subDict.getOrDefault<word>(
"fractionBasis",
"mass")
135 if (fractionBasisType ==
"mass")
140 else if (fractionBasisType ==
"mole")
147 comp[i] *= thermo_.W(i);
155 <<
"The given fractionBasis type is invalid" <<
nl 156 <<
"Valid fractionBasis types are" <<
nl 157 <<
" \"mass\" (default)" <<
nl 168 Foam::scalar Foam::functionObjects::BilgerMixtureFraction::o2Required
177 comp[i]/thermo_.W(i)*(nAtomsC_[i] + nAtomsS_[i] + 0.25*nAtomsH_[i]);
184 Foam::scalar Foam::functionObjects::BilgerMixtureFraction::o2Present
192 o2pres += comp[i]/thermo_.W(i)*nAtomsO_[i];
209 phaseName_(
dict.getOrDefault<
word>(
"phase",
word::null)),
215 IOobject::groupName(
"f_Bilger", phaseName_)
225 nSpecies_(thermo_.
Y().size()),
227 o2RequiredFuelOx_(0),
228 nAtomsC_(nSpecies_, 0),
229 nAtomsS_(nSpecies_, 0),
230 nAtomsH_(nSpecies_, 0),
231 nAtomsO_(nSpecies_, 0),
232 Yoxidiser_(nSpecies_, 0),
237 calcBilgerMixtureFraction();
260 nSpecies_ = thermo_.Y().size();
265 <<
"Number of input species is zero" 269 nAtomsC_.resize(nSpecies_, 0);
270 nAtomsS_.resize(nSpecies_, 0);
271 nAtomsH_.resize(nSpecies_, 0);
272 nAtomsO_.resize(nSpecies_, 0);
274 auto&
Y = thermo_.Y();
277 typedef BasicChemistryModel<psiReactionThermo> psiChemistryModelType;
278 typedef BasicChemistryModel<rhoReactionThermo> rhoChemistryModelType;
280 const auto* psiChemPtr =
281 mesh_.cfindObject<psiChemistryModelType>(
"chemistryProperties");
283 const auto* rhoChemPtr =
284 mesh_.cfindObject<rhoChemistryModelType>(
"chemistryProperties");
286 autoPtr<speciesCompositionTable> speciesCompPtr;
290 speciesCompPtr.reset((*psiChemPtr).thermo().specieComposition());
294 speciesCompPtr.reset((*rhoChemPtr).thermo().specieComposition());
299 <<
"BasicChemistryModel not found" 305 const List<specieElement>& curSpecieComposition =
306 (speciesCompPtr.ref())[speciesTab[i]];
308 forAll(curSpecieComposition, j)
310 const word&
e = curSpecieComposition[j].
name();
311 const label nAtoms = curSpecieComposition[j].nAtoms();
315 nAtomsC_[i] = nAtoms;
319 nAtomsS_[i] = nAtoms;
323 nAtomsH_[i] = nAtoms;
327 nAtomsO_[i] = nAtoms;
332 if (
sum(nAtomsO_) == 0)
335 <<
"No specie contains oxygen" 339 Yoxidiser_.resize(nSpecies_, 0);
340 Yfuel_.resize(nSpecies_, 0);
344 !readComposition(
dict.subDict(
"oxidiser"), Yoxidiser_)
345 || !readComposition(
dict.subDict(
"fuel"), Yfuel_)
351 o2RequiredOx_ = o2Required(Yoxidiser_) - o2Present(Yoxidiser_);
353 if (o2RequiredOx_ > 0)
356 <<
"Oxidiser composition contains not enough oxygen" <<
endl 357 <<
"Mixed up fuel and oxidiser compositions?" 361 const scalar o2RequiredFuel = o2Required(Yfuel_) - o2Present(Yfuel_);
363 if (o2RequiredFuel < 0)
366 <<
"Fuel composition contains too much oxygen" <<
endl 367 <<
"Mixed up fuel and oxidiser compositions?" 371 o2RequiredFuelOx_ = o2RequiredFuel - o2RequiredOx_;
373 if (
mag(o2RequiredFuelOx_) < SMALL)
376 <<
"Fuel and oxidiser have the same composition" 386 calcBilgerMixtureFraction();
394 return clearObject(resultName_);
401 <<
" writing field " << resultName_ <<
endl;
403 return writeObject(resultName_);
Type definitions for thermo-physics models.
defineTypeNameAndDebug(ObukhovLength, 0)
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &f1)
Abstract base-class for fluid and solid thermodynamic properties.
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
BilgerMixtureFraction(const word &name, const Time &runTime, const dictionary &dict)
Construct from Time and dictionary.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
constexpr char nl
The newline '\n' character (0x0a)
const word dictName("faMeshDefinition")
Ostream & endl(Ostream &os)
Add newline and flush stream.
Ignore writing from objectRegistry::writeObject()
const dimensionSet dimless
Dimensionless.
const Time & time() const
Return the top-level database.
Specialization of basicMultiComponentMixture for a mixture consisting of a number for molecular speci...
hashedWordList speciesTable
A table of species as a hashedWordList.
PtrList< volScalarField > & Y()
Return the mass-fraction fields.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Macros for easy insertion into run-time selection tables.
#define forAll(list, i)
Loop across all elements in list.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
const dimensionedScalar e
Elementary charge.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
A class for handling words, derived from Foam::string.
Type * getObjectPtr(const word &name, const bool recursive=false) const
Return non-const pointer to the object of the given Type, using a const-cast to have it behave like a...
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static const word null
An empty word.
virtual bool clear()
Clear the Bilger mixture-fraction field from registry.
virtual bool write()
Write Bilger mixture-fraction field.
static word timeName(const scalar t, const int precision=precision_)
Return time name of given scalar time formatted with the given precision.
static constexpr scalarRange ge0() noexcept
A greater-equals zero bound.
addToRunTimeSelectionTable(functionObject, ObukhovLength, dictionary)
const word & name() const noexcept
Return const reference to name.
virtual bool read(const dictionary &)
Read the BilgerMixtureFraction data.
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
PtrList< volScalarField > & Y
virtual bool execute()
Calculate the Bilger mixture-fraction field.
messageStream Info
Information stream (stdout output on master, null elsewhere)
virtual bool read(const dictionary &dict)
Read optional controls.
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
virtual scalar W(const label speciei) const =0
Molecular weight of the given specie [kg/kmol].
Defines the attributes of an object for which implicit objectRegistry management is supported...
Request registration (bool: true)
const fvMesh & mesh_
Reference to the fvMesh.
IOerror FatalIOError
Error stream (stdout output on all processes), with additional 'FOAM FATAL IO ERROR' header text and ...
static constexpr const zero Zero
Global zero (0)