50 void Foam::radiation::fvDOM::initialise()
53 if (mesh_.nSolutionD() == 3)
55 nRay_ = 4*nPhi_*nTheta_;
57 const scalar deltaPhi =
pi/(2.0*nPhi_);
58 const scalar deltaTheta =
pi/nTheta_;
62 for (
label m = 1; m <= 4*nPhi_; m++)
64 const scalar thetai = (2*
n - 1)*deltaTheta/2.0;
65 const scalar phii = (2*m - 1)*deltaPhi/2.0;
69 new radiativeIntensityRay
88 else if (mesh_.nSolutionD() == 2)
91 const scalar deltaTheta =
pi;
94 const scalar deltaPhi =
pi/(2.0*nPhi_);
96 for (
label m = 1; m <= 4*nPhi_; m++)
98 const scalar phii = (2*m - 1)*deltaPhi/2.0;
102 new radiativeIntensityRay
123 const scalar deltaTheta =
pi;
125 IRay_.setSize(nRay_);
126 const scalar deltaPhi =
pi;
128 for (
label m = 1; m <= 2; m++)
130 const scalar phii = (2*m - 1)*deltaPhi/2.0;
134 new radiativeIntensityRay
164 mesh_.time().timeName(),
178 if (omegaMax_ < IRay_[rayId].omega())
180 omegaMax_ = IRay_[rayId].omega();
185 Info<< typeName <<
": Created " << IRay_.size() <<
" rays with average " 186 <<
"directions (dAve) and solid angles (omega)" <<
endl;
191 <<
"Ray " << IRay_[rayId].I().
name() <<
": " 192 <<
"dAve = " << IRay_[rayId].dAve() <<
", " 193 <<
"omega = " << IRay_[rayId].omega() <<
endl;
272 nLambda_(absorptionEmission_->nBands()),
274 blackBody_(nLambda_, T),
278 coeffs_.
found(
"convergence")
280 : coeffs_.lookupOrDefault<scalar>(
"tolerance", 0)
282 maxIter_(coeffs_.lookupOrDefault<
label>(
"maxIter", 50)),
289 Foam::radiation::fvDOM::fvDOM
366 blackBody_(nLambda_, T),
413 updateBlackBodyEmission();
418 scalar maxResidual = 0;
422 Info<<
"Radiation solver iter: " << radIter <<
endl;
428 if (!rayIdConv[rayI])
430 scalar maxBandResidual = IRay_[rayI].correct();
431 maxResidual =
max(maxBandResidual, maxResidual);
433 if (maxBandResidual < tolerance_)
435 rayIdConv[rayI] =
true;
440 }
while (maxResidual > tolerance_ && radIter < maxIter_);
474 for (
label j=1; j < nLambda_; j++)
513 for (
label j=0; j < nLambda_; j++)
518 IRay_[0].ILambda(j)()*IRay_[0].omega()
521 for (
label rayI=1; rayI < nRay_; rayI++)
523 Gj.
ref() += IRay_[rayI].ILambda(j)()*IRay_[rayI].omega();
534 void Foam::radiation::fvDOM::updateBlackBodyEmission()
536 for (
label j=0; j < nLambda_; j++)
552 IRay_[rayI].addIntensity();
553 G_ += IRay_[rayI].I()*IRay_[rayI].omega();
569 const size_type i1 = name.find_first_of(
"_");
570 const size_type i2 = name.find_last_of(
"_");
dictionary coeffs_
Radiation model dictionary.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Ostream & indent(Ostream &os)
Indent stream.
A list of keyword definitions, which are a keyword followed by any number of values (e...
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
T & ref() const
Return non-const reference or generate a fatal error.
const volScalarField & T_
Reference to the temperature field.
Ostream & endl(Ostream &os)
Add newline and flush stream.
void correct(const label lambdaI, const Vector2D< scalar > &band)
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
const Time & time() const
Return the top-level database.
tmp< Foam::volScalarField > deltaLambdaT(const volScalarField &T, const Vector2D< scalar > &band) const
Proportion of total energy at T from lambda1 to lambda2.
Macros for easy insertion into run-time selection tables.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Dimension set for the base types.
void setRayIdLambdaId(const word &name, label &rayId, label &lambdaId) const
Set the rayId and lambdaId from by decomposing an intensity.
stressControl lookup("compactNormalStress") >> compactNormalStress
virtual bool read()=0
Read radiationProperties dictionary.
A class for handling words, derived from string.
void updateG()
Update G and calculate total heat flux on boundary.
virtual const fileName & name() const
Return the name of the stream.
const fvMesh & mesh_
Reference to the mesh database.
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
Top level model for radiation modelling.
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
label readLabel(Istream &is)
bool read()
Read radiation properties dictionary.
defineTypeNameAndDebug(combustionModel, 0)
Ostream & decrIndent(Ostream &os)
Decrement the indent level.
word name(const complex &)
Return a string representation of a complex.
dimensionedScalar pow3(const dimensionedScalar &ds)
virtual tmp< volScalarField::Internal > Ru() const
Source term component (constant)
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
label size_type
The type that can represent the size of a DLList.
virtual ~fvDOM()
Destructor.
Input from memory buffer stream.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
void calculate()
Solve radiation equation(s)
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
A class for managing temporary objects.
const scalar piByTwo(0.5 *pi)
Ostream & incrIndent(Ostream &os)
Increment the indent level.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
virtual tmp< volScalarField > Rp() const
Source term component (for power of T^4)
#define addToRadiationRunTimeSelectionTables(model)
autoPtr< absorptionEmissionModel > absorptionEmission_
Absorption/emission model.
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.