44namespace functionObjects
59void 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_;
100bool 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
170Foam::scalar Foam::functionObjects::BilgerMixtureFraction::o2Required
179 comp[i]/thermo_.W(i)*(nAtomsC_[i] + nAtomsS_[i] + 0.25*nAtomsH_[i]);
186Foam::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();
282 const auto* psiChemPtr =
283 mesh_.cfindObject<psiChemistryModelType>(
"chemistryProperties");
285 const auto* rhoChemPtr =
286 mesh_.cfindObject<rhoChemistryModelType>(
"chemistryProperties");
292 speciesCompPtr.
reset((*psiChemPtr).thermo().specieComposition());
296 speciesCompPtr.reset((*rhoChemPtr).thermo().specieComposition());
301 <<
"BasicChemistryModel not found"
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_);
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
Basic chemistry model templated on thermodynamics.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
virtual bool read()
Re-read model coefficients if they have changed.
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
static word timeName(const scalar t, const int precision=precision_)
Pointer management similar to std::unique_ptr, with some additional methods and type checking.
void reset(autoPtr< T > &&other) noexcept
Delete managed object and set to new given pointer.
PtrList< volScalarField > & Y()
Return the mass-fraction fields.
Specialization of basicMultiComponentMixture for a mixture consisting of a number for molecular speci...
virtual scalar W(const label speciei) const =0
Molecular weight of the given specie [kg/kmol].
Abstract base-class for fluid and solid thermodynamic properties.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Abstract base-class for Time/database function objects.
Calculates the Bilger mixture-fraction field (i.e. ) based on the elemental composition of the mixtur...
virtual bool clear()
Clear the Bilger mixture-fraction field from registry.
virtual bool execute()
Calculate the Bilger mixture-fraction field.
virtual bool write()
Write Bilger mixture-fraction field.
virtual bool read(const dictionary &)
Read the BilgerMixtureFraction data.
Specialization of Foam::functionObject for an Foam::fvMesh, providing a reference to the Foam::fvMesh...
const fvMesh & mesh_
Reference to the fvMesh.
Computes the magnitude of an input field.
const Time & time() const
Return the top-level database.
A wordList with hashed named lookup, which can be faster in some situations than using the normal lis...
bool foundObject(const word &name, const bool recursive=false) const
Is the named Type found?
Type & lookupObjectRef(const word &name, const bool recursive=false) const
static constexpr scalarRange ge0() noexcept
A greater-equals zero bound.
A class for handling words, derived from Foam::string.
static const word null
An empty word.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
PtrList< volScalarField > & Y
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const word dictName("faMeshDefinition")
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimless
Dimensionless.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
messageStream Info
Information stream (stdout output on master, null elsewhere)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
hashedWordList speciesTable
A table of species as a hashedWordList.
Ostream & endl(Ostream &os)
Add newline and flush stream.
static constexpr const zero Zero
Global zero (0)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
errorManipArg< error, int > exit(error &err, const int errNo=1)
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.
Type definitions for thermo-physics models.