39 absorptionEmissionModel,
40 wideBandAbsorptionEmission,
56 coeffsDict_((dict.
subDict(typeName +
"Coeffs"))),
58 specieIndex_(
label(0)),
61 fileName(coeffsDict_.lookup(
"lookUpTableFileName")),
80 dict.
lookup(
"bandLimits") >> iBands_[nBand];
81 dict.
lookup(
"EhrrCoeff") >> iEhrrCoeffs_[nBand];
82 totalWaveLength_ += iBands_[nBand][1] - iBands_[nBand][0];
88 const word& key = iter().keyword();
91 speciesNames_.insert(key, nSpec);
95 if (!speciesNames_.found(key))
99 "Foam::radiation::wideBandAbsorptionEmission(const" 100 "dictionary& dict, const fvMesh& mesh)" 101 ) <<
"specie: " << key <<
"is not in all the bands" 105 coeffs_[nBand][nSpec].initialise(specDicts.
subDict(key));
118 if (lookUpTable_.found(iter.key()))
120 label index = lookUpTable_.findFieldIndex(iter.key());
121 Info<<
"specie: " << iter.key() <<
" found in look-up table " 122 <<
" with index: " << index <<
endl;
123 specieIndex_[iter()] = index;
132 specieIndex_[iter()] = 0.0;
134 Info<<
"species: " << iter.key() <<
" is being solved" <<
endl;
140 "radiation::wideBandAbsorptionEmission(const" 141 "dictionary& dict, const fvMesh& mesh)" 142 ) <<
"specie: " << iter.key()
143 <<
" is neither in look-up table : " 144 << lookUpTable_.tableName() <<
" nor is being solved" 190 const List<scalar>& species = lookUpTable_.lookUp(ft[i]);
196 if (specieIndex_[
n] != 0)
199 Yipi = species[specieIndex_[
n]]*p[i]*9.869231e-6;
211 coeffs_[
n][bandI].coeffs(T[i]);
213 if (coeffs_[
n][bandI].invTemp())
221 ((((b[5]*Ti + b[4])*Ti + b[3])*Ti + b[2])*Ti + b[1])*Ti
258 if (
mesh().foundObject<volScalarField>(
"dQ"))
264 E().internalField() =
267 *(iBands_[bandI][1] - iBands_[bandI][0])
273 E().internalField() =
276 *(iBands_[bandI][1] - iBands_[bandI][0])
283 "tmp<volScalarField>" 284 "radiation::wideBandAbsorptionEmission::ECont" 289 <<
"Incompatible dimensions for dQ field" <<
endl;
305 for (
label j=0; j<nBands_; j++)
307 aLambda[j].internalField() = this->a(j);
310 aLambda[j].internalField()
311 *(iBands_[j][1] - iBands_[j][0])
dimensionedScalar pow3(const dimensionedScalar &ds)
Mesh data needed to do the Finite Volume discretisation.
const dimensionedScalar b
Wien displacement law constant: default SI units: [m.K].
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
bool foundObject(const word &name) const
Is the named Type found?
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
defineTypeNameAndDebug(cloudAbsorptionEmission, 0)
A class for handling words, derived from string.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
void size(const label)
Override size to be inconsistent with allocated storage.
errorManipArg< error, int > exit(error &err, const int errNo=1)
InternalField & internalField()
Return internal field.
const dimensionSet dimEnergy
addToRunTimeSelectionTable(absorptionEmissionModel, cloudAbsorptionEmission, dictionary)
A list of keyword definitions, which are a keyword followed by any number of values (e...
const dictionary & subDict(const word &) const
Find and return a sub-dictionary.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
const word & constant() const
Return constant name.
const Time & time() const
Return the top-level database.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Ostream & endl(Ostream &os)
Add newline and flush stream.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
#define WarningIn(functionName)
Report a warning using Foam::Warning.
Model to supply absorption and emission coefficients for radiation modelling.
PtrList< volScalarField > & Y
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Macros for easy insertion into run-time selection tables.
const dimensionSet & dimensions() const
Return dimensions.
Fundamental fluid thermodynamic properties.
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
tmp< volScalarField > aCont(const label bandI=0) const
Absorption coefficient for continuous phase.
const dimensionSet dimVolume(pow3(dimLength))
A class for handling file names.
tmp< volScalarField > ECont(const label bandI=0) const
Emission contribution for continuous phase.
virtual ~wideBandAbsorptionEmission()
Destructor.
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
wideBandAbsorptionEmission(const dictionary &dict, const fvMesh &mesh)
Construct from components.
const word dictName() const
Return the local dictionary name (final part of scoped name)
tmp< volScalarField > eCont(const label bandI=0) const
Emission coefficient for continuous phase.
void correct(volScalarField &a_, PtrList< volScalarField > &aLambda) const
Correct rays.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
A class for managing temporary objects.