Go to the documentation of this file.
51 mixedFvPatchScalarField(
p, iF),
56 valueFraction() = 1.0;
69 mixedFvPatchScalarField(ptf,
p, iF, mapper),
82 mixedFvPatchScalarField(
p, iF),
83 TName_(
dict.getOrDefault<
word>(
"T",
"T"))
85 if (
dict.found(
"refValue"))
99 valueFraction() = 1.0;
113 mixedFvPatchScalarField(ptf),
125 mixedFvPatchScalarField(ptf, iF),
148 const fvDOM& dom = db().lookupObject<
fvDOM>(
"radiationProperties");
154 const label patchi =
patch().index();
156 if (dom.nLambda() != 1)
159 <<
" a grey boundary condition is used with a non-grey "
189 const scalarField& transmissivity = ttransmissivity();
194 const vector& myRayId = dom.IRay(rayId).d();
199 for (label rayi=0; rayi < dom.nRay(); rayi++)
201 const vector& d = dom.IRay(rayi).d();
203 if ((-
n[facei] & d) < 0.0)
207 dom.IRay(rayi).ILambda(lambdaId).boundaryField()[patchi];
209 const vector& rayDave = dom.IRay(rayi).dAve();
210 Ir[facei] += IFace[facei]*(
n[facei] & rayDave);
215 if (dom.useSolarLoad())
220 dom.primaryFluxName_ +
"_0"
223 word qSecName = dom.relfectedFluxName_ +
"_0";
236 if (dom.useExternalBeam())
238 const vector sunDir = dom.solarCalc().direction();
239 const scalar directSolarRad = dom.solarCalc().directSolarRad();
243 scalar maxSunRay = -GREAT;
246 for (label rayI=0; rayI < dom.nRay(); rayI++)
248 const vector& iD = dom.IRay(rayI).d();
249 scalar dir = sunDir & iD;
257 if (rayId == SunRayId)
262 Iexternal[faceI] = directSolarRad/
mag(dom.IRay(rayId).dAve());
269 if ((-
n[faceI] & myRayId) > 0.0)
272 refGrad()[faceI] = 0.0;
273 valueFraction()[faceI] = 1.0;
275 Iexternal[faceI]*transmissivity[faceI]
277 Ir[faceI]*(scalar(1) - emissivity[faceI])
283 qem[faceI] = refValue()[faceI]*nAve[faceI];
288 valueFraction()[faceI] = 0.0;
289 refGrad()[faceI] = 0.0;
290 refValue()[faceI] = 0.0;
293 qin[faceI] = Iw[faceI]*nAve[faceI];
300 mixedFvPatchScalarField::updateCoeffs();
310 os.writeEntryIfDifferent<
word>(
"T",
"T", TName_);
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
A class for handling words, derived from Foam::string.
A class for managing temporary objects.
static constexpr const zero Zero
Global zero (0)
Different types of constants.
volScalarField & qem()
Return non-const access to the boundary emitted heat flux.
static const boundaryRadiationProperties & New(const fvMesh &mesh, Args &&... args)
Get existing or create a new MeshObject.
greyDiffusiveRadiationMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Radiation intensity for a ray in a given direction.
Unit conversion functions.
friend Ostream & operator(Ostream &, const fvPatchField< scalar > &)
This boundary condition provides a grey-diffuse condition for radiation intensity,...
const volScalarField & qr() const
Return const access to the boundary heat flux.
#define forAll(list, i)
Loop across all elements in list.
tmp< scalarField > emissivity(const label patchI, const label bandI=0, vectorField *incomingDirection=nullptr, scalarField *T=nullptr) const
Access boundary emissivity on patch.
dimensionedScalar pow4(const dimensionedScalar &ds)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
volScalarField & qin()
Return non-const access to the boundary incident heat flux.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
Boundary radiation properties holder.
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
OBJstream os(runTime.globalPath()/outputName)
Macros for easy insertion into run-time selection tables.
errorManipArg< error, int > exit(error &err, const int errNo=1)
const vector & dAve() const
Return the average vector inside the solid angle.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
static int & msgType() noexcept
Message tag of standard messages.
constexpr scalar pi(M_PI)
virtual void write(Ostream &) const
Write.
const std::string patch
OpenFOAM patch number as a std::string.
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Foam::fvPatchFieldMapper.
tmp< scalarField > transmissivity(const label patchI, const label bandI=0, vectorField *incomingDirection=nullptr, scalarField *T=nullptr) const
Access boundary transmissivity on patch.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
virtual void operator=(const UList< scalar > &)
void setRayIdLambdaId(const word &name, label &rayId, label &lambdaId) const
Set the rayId and lambdaId from by decomposing an intensity.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Finite Volume Discrete Ordinates Method. Solves the RTE equation for n directions in a participating ...