44 namespace functionObjects
50 BilgerMixtureFraction,
59 void Foam::functionObjects::BilgerMixtureFraction::calcBilgerMixtureFraction()
76 mesh_.objectRegistry::store(tCo.ptr());
81 auto&
Y = thermo_.
Y();
83 f_Bilger = -o2RequiredOx_;
88 *(nAtomsC_[i] + nAtomsS_[i] + 0.25*nAtomsH_[i] - 0.5*nAtomsO_[i])
91 f_Bilger /= o2RequiredFuelOx_;
100 bool Foam::functionObjects::BilgerMixtureFraction::readComposition
102 const dictionary& subDict,
106 auto&
Y = thermo_.Y();
113 subDict.getCheckOrDefault<scalar>
121 if (
sum(comp) < SMALL)
124 <<
"No composition is given" <<
nl 125 <<
"Valid species are:" <<
nl 132 const word fractionBasisType
134 subDict.getOrDefault<word>(
"fractionBasis",
"mass")
137 if (fractionBasisType ==
"mass")
142 else if (fractionBasisType ==
"mole")
149 comp[i] *= thermo_.W(i);
157 <<
"The given fractionBasis type is invalid" <<
nl 158 <<
"Valid fractionBasis types are" <<
nl 159 <<
" \"mass\" (default)" <<
nl 170 Foam::scalar Foam::functionObjects::BilgerMixtureFraction::o2Required
179 comp[i]/thermo_.W(i)*(nAtomsC_[i] + nAtomsS_[i] + 0.25*nAtomsH_[i]);
186 Foam::scalar Foam::functionObjects::BilgerMixtureFraction::o2Present
194 o2pres += comp[i]/thermo_.W(i)*nAtomsO_[i];
211 phaseName_(
dict.getOrDefault<
word>(
"phase",
word::null)),
217 IOobject::groupName(
"f_Bilger", phaseName_)
227 nSpecies_(thermo_.
Y().size()),
229 o2RequiredFuelOx_(0),
230 nAtomsC_(nSpecies_, 0),
231 nAtomsS_(nSpecies_, 0),
232 nAtomsH_(nSpecies_, 0),
233 nAtomsO_(nSpecies_, 0),
234 Yoxidiser_(nSpecies_, 0),
239 calcBilgerMixtureFraction();
262 nSpecies_ = thermo_.Y().size();
267 <<
"Number of input species is zero" 271 nAtomsC_.resize(nSpecies_, 0);
272 nAtomsS_.resize(nSpecies_, 0);
273 nAtomsH_.resize(nSpecies_, 0);
274 nAtomsO_.resize(nSpecies_, 0);
276 auto&
Y = thermo_.Y();
279 typedef BasicChemistryModel<psiReactionThermo> psiChemistryModelType;
280 typedef BasicChemistryModel<rhoReactionThermo> rhoChemistryModelType;
282 const auto* psiChemPtr =
283 mesh_.cfindObject<psiChemistryModelType>(
"chemistryProperties");
285 const auto* rhoChemPtr =
286 mesh_.cfindObject<rhoChemistryModelType>(
"chemistryProperties");
288 autoPtr<speciesCompositionTable> speciesCompPtr;
292 speciesCompPtr.reset((*psiChemPtr).thermo().specieComposition());
296 speciesCompPtr.reset((*rhoChemPtr).thermo().specieComposition());
301 <<
"BasicChemistryModel not found" 307 const List<specieElement>& curSpecieComposition =
308 (speciesCompPtr.ref())[speciesTab[i]];
310 forAll(curSpecieComposition, j)
312 const word&
e = curSpecieComposition[j].
name();
313 const label nAtoms = curSpecieComposition[j].nAtoms();
317 nAtomsC_[i] = nAtoms;
321 nAtomsS_[i] = nAtoms;
325 nAtomsH_[i] = nAtoms;
329 nAtomsO_[i] = nAtoms;
334 if (
sum(nAtomsO_) == 0)
337 <<
"No specie contains oxygen" 341 Yoxidiser_.resize(nSpecies_, 0);
342 Yfuel_.resize(nSpecies_, 0);
346 !readComposition(
dict.subDict(
"oxidiser"), Yoxidiser_)
347 || !readComposition(
dict.subDict(
"fuel"), Yfuel_)
353 o2RequiredOx_ = o2Required(Yoxidiser_) - o2Present(Yoxidiser_);
355 if (o2RequiredOx_ > 0)
358 <<
"Oxidiser composition contains not enough oxygen" <<
endl 359 <<
"Mixed up fuel and oxidiser compositions?" 363 const scalar o2RequiredFuel = o2Required(Yfuel_) - o2Present(Yfuel_);
365 if (o2RequiredFuel < 0)
368 <<
"Fuel composition contains too much oxygen" <<
endl 369 <<
"Mixed up fuel and oxidiser compositions?" 373 o2RequiredFuelOx_ = o2RequiredFuel - o2RequiredOx_;
375 if (
mag(o2RequiredFuelOx_) < SMALL)
378 <<
"Fuel and oxidiser have the same composition" 388 calcBilgerMixtureFraction();
396 return clearObject(resultName_);
403 <<
" writing field " << resultName_ <<
endl;
405 return writeObject(resultName_);
Type definitions for thermo-physics models.
defineTypeNameAndDebug(ObukhovLength, 0)
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)
Type & lookupObjectRef(const word &name, const bool recursive=false) const
Lookup and return non-const reference to the object of the given Type. Fatal if not found or the wron...
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.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
#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.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
static const word null
An empty word.
static tmp< T > New(Args &&... args)
Construct tmp with forwarding arguments.
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.
bool foundObject(const word &name, const bool recursive=false) const
Is the named Type found?
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...
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)