greyDiffusiveRadiationMixedFvPatchScalarField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2018 OpenFOAM Foundation
9  Copyright (C) 2016-2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
31 #include "fvPatchFieldMapper.H"
32 #include "volFields.H"
34 
35 #include "fvDOM.H"
36 #include "constants.H"
37 #include "unitConversion.H"
38 
39 using namespace Foam::constant;
40 using namespace Foam::constant::mathematical;
41 
42 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
43 
46 (
47  const fvPatch& p,
49 )
50 :
51  mixedFvPatchScalarField(p, iF),
52  TName_("T")
53 {
54  refValue() = 0.0;
55  refGrad() = 0.0;
56  valueFraction() = 1.0;
57 }
58 
59 
62 (
64  const fvPatch& p,
66  const fvPatchFieldMapper& mapper
67 )
68 :
69  mixedFvPatchScalarField(ptf, p, iF, mapper),
70  TName_(ptf.TName_)
71 {}
72 
73 
76 (
77  const fvPatch& p,
79  const dictionary& dict
80 )
81 :
82  mixedFvPatchScalarField(p, iF),
83  TName_(dict.getOrDefault<word>("T", "T"))
84 {
85  if (dict.found("refValue"))
86  {
88  (
89  scalarField("value", dict, p.size())
90  );
91  refValue() = scalarField("refValue", dict, p.size());
92  refGrad() = scalarField("refGradient", dict, p.size());
93  valueFraction() = scalarField("valueFraction", dict, p.size());
94  }
95  else
96  {
97  refValue() = 0.0;
98  refGrad() = 0.0;
99  valueFraction() = 1.0;
100 
101  fvPatchScalarField::operator=(refValue());
102  }
103 
104 }
105 
106 
109 (
111 )
112 :
113  mixedFvPatchScalarField(ptf),
114  TName_(ptf.TName_)
115 {}
116 
117 
120 (
123 )
124 :
125  mixedFvPatchScalarField(ptf, iF),
126  TName_(ptf.TName_)
127 {}
128 
129 
130 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
131 
134 {
135  if (this->updated())
136  {
137  return;
138  }
139 
140  // Since we're inside initEvaluate/evaluate there might be processor
141  // comms underway. Change the tag we use.
142  int oldTag = UPstream::msgType();
143  UPstream::msgType() = oldTag+1;
144 
145  const scalarField& Tp =
146  patch().lookupPatchField<volScalarField, scalar>(TName_);
147 
148  const fvDOM& dom = db().lookupObject<fvDOM>("radiationProperties");
149 
150  label rayId = -1;
151  label lambdaId = -1;
152  dom.setRayIdLambdaId(internalField().name(), rayId, lambdaId);
153 
154  const label patchi = patch().index();
155 
156  if (dom.nLambda() != 1)
157  {
159  << " a grey boundary condition is used with a non-grey "
160  << "absorption model" << nl << exit(FatalError);
161  }
162 
163  scalarField& Iw = *this;
164 
165  const vectorField n(patch().nf());
166 
167  radiativeIntensityRay& ray =
168  const_cast<radiativeIntensityRay&>(dom.IRay(rayId));
169 
170  const scalarField nAve(n & ray.dAve());
171 
172  ray.qr().boundaryFieldRef()[patchi] += Iw*nAve;
173 
174  const boundaryRadiationProperties& boundaryRadiation =
175  boundaryRadiationProperties::New(internalField().mesh());
176 
177  const tmp<scalarField> temissivity
178  (
179  boundaryRadiation.emissivity(patch().index())
180  );
181 
182  const scalarField& emissivity = temissivity();
183 
184  const tmp<scalarField> ttransmissivity
185  (
186  boundaryRadiation.transmissivity(patch().index())
187  );
188 
189  const scalarField& transmissivity = ttransmissivity();
190 
191  scalarField& qem = ray.qem().boundaryFieldRef()[patchi];
192  scalarField& qin = ray.qin().boundaryFieldRef()[patchi];
193 
194  const vector& myRayId = dom.IRay(rayId).d();
195 
196  scalarField Ir(patch().size(), Zero);
197  forAll(Iw, facei)
198  {
199  for (label rayi=0; rayi < dom.nRay(); rayi++)
200  {
201  const vector& d = dom.IRay(rayi).d();
202 
203  if ((-n[facei] & d) < 0.0)
204  {
205  // q into the wall
206  const scalarField& IFace =
207  dom.IRay(rayi).ILambda(lambdaId).boundaryField()[patchi];
208 
209  const vector& rayDave = dom.IRay(rayi).dAve();
210  Ir[facei] += IFace[facei]*(n[facei] & rayDave);
211  }
212  }
213  }
214 
215  if (dom.useSolarLoad())
216  {
217  // Looking for primary heat flux single band
218  Ir += patch().lookupPatchField<volScalarField,scalar>
219  (
220  dom.primaryFluxName_ + "_0"
221  );
222 
223  word qSecName = dom.relfectedFluxName_ + "_0";
224 
225  if (this->db().foundObject<volScalarField>(qSecName))
226  {
227  const volScalarField& qSec =
228  this->db().lookupObject<volScalarField>(qSecName);
229 
230  Ir += qSec.boundaryField()[patch().index()];
231  }
232  }
233 
234  scalarField Iexternal(this->size(), 0.0);
235 
236  if (dom.useExternalBeam())
237  {
238  const vector sunDir = dom.solarCalc().direction();
239  const scalar directSolarRad = dom.solarCalc().directSolarRad();
240 
241  //label nRaysBeam = dom.nRaysBeam();
242  label SunRayId(-1);
243  scalar maxSunRay = -GREAT;
244 
245  // Looking for the ray closest to the Sun direction
246  for (label rayI=0; rayI < dom.nRay(); rayI++)
247  {
248  const vector& iD = dom.IRay(rayI).d();
249  scalar dir = sunDir & iD;
250  if (dir > maxSunRay)
251  {
252  maxSunRay = dir;
253  SunRayId = rayI;
254  }
255  }
256 
257  if (rayId == SunRayId)
258  {
259  const scalarField nAve(n & dom.IRay(rayId).dAve());
260  forAll(Iexternal, faceI)
261  {
262  Iexternal[faceI] = directSolarRad/mag(dom.IRay(rayId).dAve());
263  }
264  }
265  }
266 
267  forAll(Iw, faceI)
268  {
269  if ((-n[faceI] & myRayId) > 0.0)
270  {
271  // direction out of the wall
272  refGrad()[faceI] = 0.0;
273  valueFraction()[faceI] = 1.0;
274  refValue()[faceI] =
275  Iexternal[faceI]*transmissivity[faceI]
276  + (
277  Ir[faceI]*(scalar(1) - emissivity[faceI])
278  + emissivity[faceI]*physicoChemical::sigma.value()
279  * pow4(Tp[faceI])
280  )/pi;
281 
282  // Emitted heat flux from this ray direction
283  qem[faceI] = refValue()[faceI]*nAve[faceI];
284  }
285  else
286  {
287  // direction into the wall
288  valueFraction()[faceI] = 0.0;
289  refGrad()[faceI] = 0.0;
290  refValue()[faceI] = 0.0; //not used
291 
292  // Incident heat flux on this ray direction
293  qin[faceI] = Iw[faceI]*nAve[faceI];
294  }
295  }
296 
297  // Restore tag
298  UPstream::msgType() = oldTag;
299 
300  mixedFvPatchScalarField::updateCoeffs();
301 }
302 
303 
305 (
306  Ostream& os
307 ) const
308 {
310  os.writeEntryIfDifferent<word>("T", "T", TName_);
311 }
312 
313 
314 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
315 
316 namespace Foam
317 {
318 namespace radiation
319 {
321  (
324  );
325 }
326 }
327 
328 
329 // ************************************************************************* //
Foam::fvPatchField< scalar >
volFields.H
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:52
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::constant
Different types of constants.
Definition: atomicConstants.C:38
Foam::radiation::radiativeIntensityRay::qem
volScalarField & qem()
Return non-const access to the boundary emitted heat flux.
Definition: radiativeIntensityRayI.H:67
Foam::MeshObject< fvMesh, Foam::GeometricMeshObject, boundaryRadiationProperties >::New
static const boundaryRadiationProperties & New(const fvMesh &mesh, Args &&... args)
Get existing or create a new MeshObject.
Definition: MeshObject.C:48
Foam::radiation::greyDiffusiveRadiationMixedFvPatchScalarField::greyDiffusiveRadiationMixedFvPatchScalarField
greyDiffusiveRadiationMixedFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Definition: greyDiffusiveRadiationMixedFvPatchScalarField.C:46
Foam::radiation::radiativeIntensityRay
Radiation intensity for a ray in a given direction.
Definition: radiativeIntensityRay.H:57
unitConversion.H
Unit conversion functions.
Foam::fvPatchField< scalar >::operator
friend Ostream & operator(Ostream &, const fvPatchField< scalar > &)
fvPatchFieldMapper.H
Foam::radiation::greyDiffusiveRadiationMixedFvPatchScalarField
This boundary condition provides a grey-diffuse condition for radiation intensity,...
Definition: greyDiffusiveRadiationMixedFvPatchScalarField.H:89
Foam::radiation::radiativeIntensityRay::qr
const volScalarField & qr() const
Return const access to the boundary heat flux.
Definition: radiativeIntensityRayI.H:36
fvDOM.H
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::radiation::boundaryRadiationProperties::emissivity
tmp< scalarField > emissivity(const label patchI, const label bandI=0, vectorField *incomingDirection=nullptr, scalarField *T=nullptr) const
Access boundary emissivity on patch.
Definition: boundaryRadiationProperties.C:111
n
label n
Definition: TABSMDCalcMethod2.H:31
Foam::pow4
dimensionedScalar pow4(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:100
Foam::Field< scalar >
Foam::fvPatch
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:65
Foam::radiation::radiativeIntensityRay::qin
volScalarField & qin()
Return non-const access to the boundary incident heat flux.
Definition: radiativeIntensityRayI.H:54
Foam::volScalarField
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
radiation
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
Foam::radiation::boundaryRadiationProperties
Boundary radiation properties holder.
Definition: boundaryRadiationProperties.H:56
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::FatalError
error FatalError
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
os
OBJstream os(runTime.globalPath()/outputName)
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
constants.H
Foam::exit
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
Foam::radiation::radiativeIntensityRay::dAve
const vector & dAve() const
Return the average vector inside the solid angle.
Definition: radiativeIntensityRayI.H:79
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::constant::mathematical
Mathematical constants.
Foam::UPstream::msgType
static int & msgType() noexcept
Message tag of standard messages.
Definition: UPstream.H:540
Foam::constant::mathematical::pi
constexpr scalar pi(M_PI)
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::radiation::greyDiffusiveRadiationMixedFvPatchScalarField::write
virtual void write(Ostream &) const
Write.
Definition: greyDiffusiveRadiationMixedFvPatchScalarField.C:305
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
Foam::GeometricField::boundaryFieldRef
Boundary & boundaryFieldRef(const bool updateAccessTime=true)
Return a reference to the boundary field.
Definition: GeometricField.C:783
Foam::Vector< scalar >
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::radiation::greyDiffusiveRadiationMixedFvPatchScalarField::updateCoeffs
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: greyDiffusiveRadiationMixedFvPatchScalarField.C:133
Foam::constant::physicoChemical::sigma
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
Foam::vtk::write
void write(vtk::formatter &fmt, const Type &val, const label n=1)
Component-wise write of a value (N times)
Definition: foamVtkOutputTemplates.C:36
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
Foam::fvPatchFieldMapper
Foam::fvPatchFieldMapper.
Definition: fvPatchFieldMapper.H:47
Foam::radiation::boundaryRadiationProperties::transmissivity
tmp< scalarField > transmissivity(const label patchI, const label bandI=0, vectorField *incomingDirection=nullptr, scalarField *T=nullptr) const
Access boundary transmissivity on patch.
Definition: boundaryRadiationProperties.C:229
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::makePatchTypeField
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
greyDiffusiveRadiationMixedFvPatchScalarField.H
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::fvPatchField< scalar >::operator=
virtual void operator=(const UList< scalar > &)
Definition: fvPatchField.C:404
Foam::radiation::fvDOM::setRayIdLambdaId
void setRayIdLambdaId(const word &name, label &rayId, label &lambdaId) const
Set the rayId and lambdaId from by decomposing an intensity.
Definition: fvDOM.C:761
boundaryRadiationProperties.H
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62
Foam::DimensionedField
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Definition: DimensionedField.H:54
Foam::radiation::fvDOM
Finite Volume Discrete Ordinates Method. Solves the RTE equation for n directions in a participating ...
Definition: fvDOM.H:118