51 mixedFvPatchScalarField(
p, iF),
58 valueFraction() = 1.0;
71 mixedFvPatchScalarField(ptf,
p, iF, mapper),
73 qRadExt_(ptf.qRadExt_),
74 qRadExtDir_(ptf.qRadExtDir_)
86 mixedFvPatchScalarField(
p, iF),
87 TName_(
dict.getOrDefault<
word>(
"T",
"T")),
88 qRadExt_(
dict.getOrDefault<scalar>(
"qRadExt", 0)),
91 if (
dict.found(
"refValue"))
105 valueFraction() = 1.0;
118 mixedFvPatchScalarField(ptf),
120 qRadExt_(ptf.qRadExt_),
121 qRadExtDir_(ptf.qRadExtDir_)
132 mixedFvPatchScalarField(ptf, iF),
134 qRadExt_(ptf.qRadExt_),
135 qRadExtDir_(ptf.qRadExtDir_)
157 const fvDOM& dom = db().lookupObject<
fvDOM>(
"radiationProperties");
163 const label patchi =
patch().index();
165 if (dom.nLambda() != 1)
168 <<
" a grey boundary condition is used with a non-grey " 176 radiativeIntensityRay& ray =
177 const_cast<radiativeIntensityRay&
>(dom.IRay(rayId));
181 ray.qr().boundaryFieldRef()[patchi] += Iw*nAve;
183 const boundaryRadiationProperties& boundaryRadiation =
186 const tmp<scalarField> temissivity
188 boundaryRadiation.emissivity(
patch().index())
193 const tmp<scalarField> ttransmissivity
195 boundaryRadiation.transmissivity(
patch().index())
198 const scalarField& transmissivity = ttransmissivity();
200 scalarField& qem = ray.qem().boundaryFieldRef()[patchi];
201 scalarField& qin = ray.qin().boundaryFieldRef()[patchi];
203 const vector& myRayId = dom.IRay(rayId).d();
208 for (label rayi=0; rayi < dom.nRay(); rayi++)
210 const vector& d = dom.IRay(rayi).d();
212 if ((-
n[facei] & d) < 0.0)
216 dom.IRay(rayi).ILambda(lambdaId).boundaryField()[patchi];
218 const vector& rayDave = dom.IRay(rayi).dAve();
219 Ir[facei] += IFace[facei]*(
n[facei] & rayDave);
224 if (dom.useSolarLoad())
229 dom.primaryFluxName_ +
"_0" 232 word qSecName = dom.relfectedFluxName_ +
"_0";
245 if (dom.useExternalBeam())
247 const vector sunDir = dom.solarCalc().direction();
248 const scalar directSolarRad = dom.solarCalc().directSolarRad();
252 scalar maxSunRay = -GREAT;
255 for (label rayI=0; rayI < dom.nRay(); rayI++)
257 const vector& iD = dom.IRay(rayI).d();
258 scalar dir = sunDir & iD;
266 if (rayId == SunRayId)
271 Iexternal[faceI] = directSolarRad/
mag(dom.IRay(rayId).dAve());
280 if (
mag(qRadExtDir_) > 0)
283 scalar maxRay = -GREAT;
286 for (label rayI = 0; rayI < dom.nRay(); ++rayI)
288 const vector& iD = dom.IRay(rayI).d();
289 const scalar dir = qRadExtDir_ & iD;
298 if (rayId == rayqoId)
302 Isource[faceI] += qRadExt_/
mag(dom.IRay(rayId).dAve());
311 scalar maxRay = -GREAT;
314 for (label rayI = 0; rayI < dom.nRay(); ++rayI)
316 const vector& iD = dom.IRay(rayI).d();
317 const scalar dir = -
n[faceI] & iD;
326 if (rayId == rayqoId)
328 Isource[faceI] += qRadExt_/
mag(dom.IRay(rayId).dAve());
336 if ((-
n[faceI] & myRayId) > 0.0)
339 refGrad()[faceI] = 0.0;
340 valueFraction()[faceI] = 1.0;
343 + Iexternal[faceI]*transmissivity[faceI]
345 Ir[faceI]*(scalar(1) - emissivity[faceI])
351 qem[faceI] = refValue()[faceI]*nAve[faceI];
356 valueFraction()[faceI] = 0.0;
357 refGrad()[faceI] = 0.0;
358 refValue()[faceI] = 0.0;
361 qin[faceI] = Iw[faceI]*nAve[faceI];
368 mixedFvPatchScalarField::updateCoeffs();
378 os.writeEntryIfDifferent<
word>(
"T",
"T", TName_);
379 os.writeEntryIfDifferent<scalar>(
"qRadExt",
Zero, qRadExt_);
380 os.writeEntryIfDifferent<
vector>(
"qRadExtDir",
Zero, qRadExtDir_);
Different types of constants.
virtual void write(Ostream &) const
Write.
makePatchTypeField(fvPatchScalarField, greyDiffusiveRadiationMixedFvPatchScalarField)
errorManipArg< error, int > exit(error &err, const int errNo=1)
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
error FatalError
Error stream (stdout output on all processes), with additional 'FOAM FATAL ERROR' header text and sta...
A list of keyword definitions, which are a keyword followed by a number of values (eg...
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
This boundary condition provides a grey-diffuse condition for radiation intensity, I, for use with the finite-volume discrete-ordinates model (fvDOM), in which the radiation temperature is retrieved from the temperature field boundary condition.
friend Ostream & operator(Ostream &, const fvPatchField< scalar > &)
static const boundaryRadiationProperties & New(const fvMesh &mesh, Args &&... args)
Get existing or create a new MeshObject.
Unit conversion functions.
constexpr char nl
The newline '\n' character (0x0a)
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
static int & msgType() noexcept
Message tag of standard messages.
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
Macros for easy insertion into run-time selection tables.
#define forAll(list, i)
Loop across all elements in list.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
void setRayIdLambdaId(const word &name, label &rayId, label &lambdaId) const
Set the rayId and lambdaId from by decomposing an intensity.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
A class for handling words, derived from Foam::string.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
A FieldMapper for finite-volume patch fields.
constexpr scalar pi(M_PI)
greyDiffusiveRadiationMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
OBJstream os(runTime.globalPath()/outputName)
virtual void operator=(const UList< scalar > &)
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensionedScalar pow4(const dimensionedScalar &ds)
const std::string patch
OpenFOAM patch number as a std::string.
Field< vector > vectorField
Specialisation of Field<T> for vector.
Finite Volume Discrete Ordinates Method. Solves the RTE equation for n directions in a participating ...
const Boundary & boundaryField() const noexcept
Return const-reference to the boundary field.
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
static constexpr const zero Zero
Global zero (0)