greyDiffusiveRadiationMixedFvPatchScalarField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2024 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
28 #include "fieldMapper.H"
29 #include "volFields.H"
30 
31 #include "fvDOM.H"
32 #include "constants.H"
33 
34 using namespace Foam::constant;
35 using namespace Foam::constant::mathematical;
36 
37 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
38 
41 (
42  const fvPatch& p,
44  const dictionary& dict
45 )
46 :
47  mixedFvPatchScalarField(p, iF, dict, false),
49  TName_(dict.lookupOrDefault<word>("T", "T"))
50 {
51  if (dict.found("refValue"))
52  {
53  refValue() = scalarField("refValue", iF.dimensions(), dict, p.size());
54  refGrad() =
56  (
57  "refGradient",
58  iF.dimensions()/dimLength,
59  dict,
60  p.size()
61  );
62  valueFraction() =
63  scalarField("valueFraction", unitFraction, dict, p.size());
64 
66  (
67  scalarField("value", iF.dimensions(), dict, p.size())
68  );
69  }
70  else
71  {
72  refValue() = 0.0;
73  refGrad() = 0.0;
74  valueFraction() = 1.0;
75 
77  }
78 }
79 
80 
83 (
85  const fvPatch& p,
87  const fieldMapper& mapper
88 )
89 :
90  mixedFvPatchScalarField(ptf, p, iF, mapper),
92  (
93  p,
94  ptf.emissivityMethod(),
95  ptf.emissivity_,
96  mapper
97  ),
98  TName_(ptf.TName_)
99 {}
100 
101 
104 (
107 )
108 :
109  mixedFvPatchScalarField(ptf, iF),
111  (
112  ptf.patch(),
113  ptf.emissivityMethod(),
114  ptf.emissivity_
115  ),
116  TName_(ptf.TName_)
117 {}
118 
119 
120 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
121 
123 (
124  const fvPatchScalarField& ptf,
125  const fieldMapper& mapper
126 )
127 {
128  mixedFvPatchScalarField::map(ptf, mapper);
129  radiationCoupledBase::map(ptf, mapper);
130 }
131 
132 
134 (
135  const fvPatchScalarField& ptf
136 )
137 {
138  mixedFvPatchScalarField::reset(ptf);
140 }
141 
142 
144 {
145  if (this->updated())
146  {
147  return;
148  }
149 
150  // Since we're inside initEvaluate/evaluate there might be processor
151  // comms underway. Change the tag we use.
152  int oldTag = UPstream::msgType();
153  UPstream::msgType() = oldTag+1;
154 
155  const scalarField& Tp =
156  patch().lookupPatchField<volScalarField, scalar>(TName_);
157 
158  const radiationModels::fvDOM& dom =
159  db().lookupObject<radiationModels::fvDOM>("radiationProperties");
160 
161  label rayId = -1;
162  label lambdaId = -1;
163  dom.setRayIdLambdaId(internalField().name(), rayId, lambdaId);
164 
165  const label patchi = patch().index();
166 
167  if (dom.nLambda() != 1)
168  {
170  << " a grey boundary condition is used with a non-grey "
171  << "absorption model" << nl << exit(FatalError);
172  }
173 
174  scalarField& Iw = *this;
175  const vectorField n(patch().nf());
176 
178  const_cast<radiationModels::radiativeIntensityRay&>(dom.IRay(rayId));
179 
180  const scalarField nAve(n & ray.dAve());
181 
182  ray.qr().boundaryFieldRef()[patchi] += Iw*nAve;
183 
184  const scalarField emissivity(this->emissivity());
185 
186  scalarField& qem = ray.qem().boundaryFieldRef()[patchi];
187  scalarField& qin = ray.qin().boundaryFieldRef()[patchi];
188 
189  const vector& myRayId = dom.IRay(rayId).dAve();
190 
191  // Use updated Ir while iterating over rays
192  // avoids to used lagged qin
193  scalarField Ir = dom.IRay(0).qin().boundaryField()[patchi];
194 
195  for (label rayI=1; rayI < dom.nRay(); rayI++)
196  {
197  Ir += dom.IRay(rayI).qin().boundaryField()[patchi];
198  }
199 
200  forAll(Iw, facei)
201  {
202  if ((-n[facei] & myRayId) > 0.0)
203  {
204  // direction out of the wall
205  refGrad()[facei] = 0.0;
206  valueFraction()[facei] = 1.0;
207  refValue()[facei] =
208  (
209  Ir[facei]*(scalar(1) - emissivity[facei])
210  + emissivity[facei]*physicoChemical::sigma.value()
211  * pow4(Tp[facei])
212  )/pi;
213 
214  // Emitted heat flux from this ray direction
215  qem[facei] = refValue()[facei]*nAve[facei];
216  }
217  else
218  {
219  // direction into the wall
220  valueFraction()[facei] = 0.0;
221  refGrad()[facei] = 0.0;
222  refValue()[facei] = 0.0; // not used
223 
224  // Incident heat flux on this ray direction
225  qin[facei] = Iw[facei]*nAve[facei];
226  }
227  }
228 
229  // Restore tag
230  UPstream::msgType() = oldTag;
231 
232  mixedFvPatchScalarField::updateCoeffs();
233 }
234 
235 
237 (
238  Ostream& os
239 ) const
240 {
243  writeEntryIfDifferent<word>(os, "T", "T", TName_);
244 }
245 
246 
247 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
248 
249 namespace Foam
250 {
252  (
255  );
256 }
257 
258 
259 // ************************************************************************* //
label n
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet & dimensions() const
Return dimensions.
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:476
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Abstract base class for field mapping.
Definition: fieldMapper.H:48
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:88
friend Ostream & operator(Ostream &, const fvPatchField< Type > &)
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:249
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
This boundary condition provides a grey-diffuse condition for radiation intensity,...
greyDiffusiveRadiationMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
virtual void reset(const fvPatchScalarField &)
Reset the fvPatchField to the given fvPatchField.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual void map(const fvPatchScalarField &, const fieldMapper &)
Map the given fvPatchField onto this fvPatchField.
Common functions to emissivity. It gets supplied from lookup into a dictionary or calculated by the s...
void write(Ostream &) const
Write.
virtual void reset(const fvPatchScalarField &)
Reset the fvPatchField to the given fvPatchField.
virtual void map(const fvPatchScalarField &, const fieldMapper &)
Map the given fvPatchField onto this fvPatchField.
Finite Volume Discrete Ordinates Method. Solves the RTE equation for n directions in a participating ...
Definition: fvDOM.H:77
void setRayIdLambdaId(const word &name, label &rayId, label &lambdaId) const
Set the rayId and lambdaId from by decomposing an intensity.
Definition: fvDOM.C:577
Radiation intensity for a ray in a given direction.
const vector & dAve() const
Return the average vector inside the solid angle.
const volScalarField & qr() const
Return const access to the boundary heat flux.
volScalarField & qin()
Return non-const access to the boundary incident heat flux.
volScalarField & qem()
Return non-const access to the boundary emitted heat flux.
A class for handling words, derived from string.
Definition: word.H:62
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:334
label patchi
mathematical constants.
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m^2/K^4].
Collection of constants.
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
Namespace for OpenFOAM.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
const dimensionSet dimLength
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
dimensionedScalar pow4(const dimensionedScalar &ds)
error FatalError
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
static const char nl
Definition: Ostream.H:266
const unitConversion unitFraction
dictionary dict
volScalarField & p