Go to the documentation of this file.
43 absorptionEmissionModel,
44 wideBandAbsorptionEmission,
60 coeffsDict_((
dict.optionalSubDict(typeName +
"Coeffs"))),
69 const dictionary& functionDicts =
dict.optionalSubDict(typeName +
"Coeffs");
70 for (
const entry& dEntry : functionDicts)
79 dict.readEntry(
"bandLimits", iBands_[nBand]);
80 dict.readEntry(
"EhrrCoeff", iEhrrCoeffs_[nBand]);
81 totalWaveLength_ += iBands_[nBand][1] - iBands_[nBand][0];
86 for (
const entry& dEntry : specDicts)
88 const word& key = dEntry.keyword();
92 speciesNames_.insert(key, nSpec);
94 else if (!speciesNames_.found(key))
97 <<
"specie: " << key <<
" is not in all the bands"
100 coeffs_[nBand][nSpec].initialise(specDicts.subDict(key));
109 coeffsDict_.found(
"lookUpTableFileName")
110 &&
"none" != coeffsDict_.get<
word>(
"lookUpTableFileName")
117 coeffsDict_.get<
fileName>(
"lookUpTableFileName"),
118 mesh.time().constant(),
126 <<
"specie ft is not present to use with "
127 <<
"lookUpTableFileName " <<
nl
138 const word& specieName = iter.key();
139 const label index = iter.val();
143 if (!lookUpTablePtr_.empty())
145 if (lookUpTablePtr_().found(specieName))
147 const label fieldIndex =
148 lookUpTablePtr_().findFieldIndex(specieName);
150 Info<<
"specie: " << specieName <<
" found on look-up table "
151 <<
" with index: " << fieldIndex <<
endl;
153 specieIndex_[index] = fieldIndex;
158 specieIndex_[index] = 0;
160 Info<<
"specie: " << specieName <<
" is being solved" <<
endl;
165 <<
"specie: " << specieName
166 <<
" is neither in look-up table: "
167 << lookUpTablePtr_().tableName()
168 <<
" nor is being solved" <<
nl
175 specieIndex_[index] = 0;
181 <<
"There is no lookup table and the specie" <<
nl
183 <<
" is not found " <<
nl
231 const label
n = iter();
233 if (specieIndex_[
n] != 0)
238 const List<scalar>& Ynft = lookUpTablePtr_().lookUp(ft[celli]);
241 Xipi = Ynft[specieIndex_[
n]]*
paToAtm(
p[celli]);
251 const label index =
mixture.species()[iter.key()];
259 scalar Ti =
T[celli];
262 coeffs_[bandi][
n].coeffs(
T[celli]);
264 if (coeffs_[bandi][
n].invTemp())
272 ((((
b[5]*Ti +
b[4])*Ti +
b[3])*Ti +
b[2])*Ti +
b[1])*Ti
319 *
Qdot.primitiveField()
320 *(iBands_[bandi][1] - iBands_[bandi][0])
328 *
Qdot.primitiveField()
329 *(iBands_[bandi][1] - iBands_[bandi][0])
335 <<
"Incompatible dimensions for Qdot field" <<
endl;
351 for (label j=0; j<nBands_; j++)
353 aLambda[j].primitiveFieldRef() = this->a(j);
356 aLambda[j].primitiveField()
357 *(iBands_[j][1] - iBands_[j][0])
A keyword and a list of tokens is an 'entry'.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Dimensionless.
A class for handling words, derived from Foam::string.
A class for handling file names.
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< " ";}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< " ";}gmvFile<< nl;forAll(lagrangianScalarNames, i){ word name=lagrangianScalarNames[i];IOField< scalar > s(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Specialization of basicMultiComponentMixture for a mixture consisting of a number for molecular speci...
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
const dimensionSet dimEnergy
tmp< volScalarField > ECont(const label bandi=0) const
Emission contribution for continuous phase.
Fundamental fluid thermodynamic properties.
const word dictName("blockMeshDict")
Unit conversion functions.
Ostream & endl(Ostream &os)
Add newline and flush stream.
#define forAll(list, i)
Loop across all elements in list.
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
dimensionedScalar pow3(const dimensionedScalar &ds)
constexpr scalar paToAtm(const scalar pa) noexcept
Conversion from Pa to atm.
messageStream Info
Information stream (uses stdout - output is on the master only)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
virtual ~wideBandAbsorptionEmission()
Destructor.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers....
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Macros for easy insertion into run-time selection tables.
Mesh data needed to do the Finite Volume discretisation.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
defineTypeNameAndDebug(cloudAbsorptionEmission, 0)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Internal::FieldType & primitiveFieldRef(const bool updateAccessTime=true)
Return a reference to the internal field.
wideBandAbsorptionEmission(const dictionary &dict, const fvMesh &mesh)
Construct from components.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
tmp< volScalarField > eCont(const label bandi=0) const
Emission coefficient for continuous phase.
forAllConstIters(mixture.phases(), phase)
tmp< volScalarField > aCont(const label bandi=0) const
Absorption coefficient for continuous phase.
Model to supply absorption and emission coefficients for radiation modelling.
const dimensionSet dimVolume(pow3(dimLength))
#define WarningInFunction
Report a warning using Foam::Warning.
addToRunTimeSelectionTable(absorptionEmissionModel, cloudAbsorptionEmission, dictionary)
void correct(volScalarField &a, PtrList< volScalarField > &aLambda) const
Correct rays.