56 const volScalarField&
T = thermo_.T();
57 const volScalarField&
p = thermo_.p();
59 tmp<scalarField> tXj(
new scalarField(
T.primitiveField().size(), Zero));
60 scalarField& Xj = tXj.ref();
62 tmp<scalarField> tRhoInv(
new scalarField(
T.primitiveField().size(), Zero));
63 scalarField& rhoInv = tRhoInv.ref();
65 forAll(mixture_.Y(), specieI)
67 const scalarField& Yi = mixture_.Y()[specieI];
72 Yi[iCell]/mixture_.rho(specieI,
p[iCell],
T[iCell]);
76 const label mySpecieI = mixture_.species().
find(specie);
79 Xj[iCell] = Yj[iCell]/mixture_.rho(mySpecieI,
p[iCell],
T[iCell]);
87Foam::radiation::greyMeanSolidAbsorptionEmission::
88greyMeanSolidAbsorptionEmission
95 coeffsDict_((
dict.optionalSubDict(typeName +
"Coeffs"))),
99 solidData_(mixture_.
Y().size())
101 if (!isA<basicSpecieMixture>(thermo_))
104 <<
"Model requires a multi-component thermo package"
111 for (
const entry& dEntry : functionDicts)
113 if (!dEntry.isDict())
118 const word& key = dEntry.keyword();
124 <<
" specie: " << key <<
" is not found in the solid mixture"
126 <<
" specie is the mixture are:" << mixture_.
species() <<
nl
129 speciesNames_.
insert(key, nFunc);
131 dict.
readEntry(
"absorptivity", solidData_[nFunc][absorptivity]);
142Foam::radiation::greyMeanSolidAbsorptionEmission::
143calc(
const label propertyId)
const
159 extrapolatedCalculatedFvPatchVectorField::typeName
167 if (mixture_.contains(iter.key()))
169 a += solidData_[iter.val()][propertyId]*X(iter.key());
173 ta.ref().correctBoundaryConditions();
184 return calc(emissivity);
194 return calc(absorptivity);
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
bool insert(const Key &key, const T &obj)
Copy insert a new entry, not overwriting existing entries.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
label find(const T &val, label pos=0) const
Find index of the first occurrence of the value.
const speciesTable & species() const
Return the table of species.
bool contains(const word &specieName) const
Does the mixture include this specie?
Specialization of basicMultiComponentMixture for a mixture consisting of a number for molecular speci...
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,...
const dictionary & optionalSubDict(const word &keyword, enum keyType::option matchOpt=keyType::REGEX) const
Find and return a sub-dictionary, otherwise return this dictionary.
bool readEntry(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX, bool mandatory=true) const
A keyword and a list of tokens is an 'entry'.
Mesh data needed to do the Finite Volume discretisation.
Model to supply absorption and emission coefficients for radiation modelling.
const dictionary & dict() const
Reference to the dictionary.
greyMeanSolidAbsorptionEmission radiation absorption and emission coefficients for solid mixture
tmp< volScalarField > eCont(const label bandI=0) const
Emission coefficient for continuous phase.
tmp< volScalarField > aCont(const label bandI=0) const
Absorption coefficient for continuous phase.
Fundamental solid thermodynamic properties.
A class for managing temporary objects.
A class for handling words, derived from Foam::string.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
PtrList< volScalarField > & Y
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
const word dictName("faMeshDefinition")
#define WarningInFunction
Report a warning using Foam::Warning.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
const dimensionSet dimless
Dimensionless.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Ostream & endl(Ostream &os)
Add newline and flush stream.
errorManip< error > abort(error &err)
static constexpr const zero Zero
Global zero (0)
constexpr char nl
The newline '\n' character (0x0a)
#define forAll(list, i)
Loop across all elements in list.
#define forAllConstIters(container, iter)
Iterate across all elements of the container object with const access.
Unit conversion functions.