MarshakRadiationFvPatchScalarField.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 #include "radiationModel.H"
32 
33 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
34 
36 (
37  const fvPatch& p,
39  const dictionary& dict
40 )
41 :
42  mixedFvPatchScalarField(p, iF, dict, false),
44  TName_(dict.lookupOrDefault<word>("T", "T"))
45 {
46  if (dict.found("value"))
47  {
48  refValue() = scalarField("value", iF.dimensions(), dict, p.size());
49  }
50  else
51  {
52  refValue() = 0.0;
53  }
54 
55  // zero gradient
56  refGrad() = 0.0;
57 
58  valueFraction() = 1.0;
59 
61 }
62 
63 
65 (
67  const fvPatch& p,
69  const fieldMapper& mapper
70 )
71 :
72  mixedFvPatchScalarField(ptf, p, iF, mapper),
74  (
75  p,
76  ptf.emissivityMethod(),
77  ptf.emissivity_,
78  mapper
79  ),
80  TName_(ptf.TName_)
81 {}
82 
83 
85 (
88 )
89 :
90  mixedFvPatchScalarField(ptf, iF),
92  (
93  ptf.patch(),
94  ptf.emissivityMethod(),
95  ptf.emissivity_
96  ),
97  TName_(ptf.TName_)
98 {}
99 
100 
101 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
102 
104 (
105  const fvPatchScalarField& ptf,
106  const fieldMapper& mapper
107 )
108 {
109  mixedFvPatchScalarField::map(ptf, mapper);
110  radiationCoupledBase::map(ptf, mapper);
111 }
112 
113 
115 (
116  const fvPatchScalarField& ptf
117 )
118 {
119  mixedFvPatchScalarField::reset(ptf);
121 }
122 
123 
125 {
126  if (this->updated())
127  {
128  return;
129  }
130 
131  // Since we're inside initEvaluate/evaluate there might be processor
132  // comms underway. Change the tag we use.
133  int oldTag = UPstream::msgType();
134  UPstream::msgType() = oldTag+1;
135 
136  // Temperature field
137  const scalarField& Tp =
138  patch().lookupPatchField<volScalarField, scalar>(TName_);
139 
140  // Re-calc reference value
141  refValue() = 4.0*constant::physicoChemical::sigma.value()*pow4(Tp);
142 
143  // Diffusion coefficient - created by radiation model's ::updateCoeffs()
144  const scalarField& gamma =
145  patch().lookupPatchField<volScalarField, scalar>("gammaRad");
146 
147  const scalarField emissivity(this->emissivity());
148 
149  const scalarField Ep(emissivity/(2.0*(2.0 - emissivity)));
150 
151  // Set value fraction
152  valueFraction() = 1.0/(1.0 + gamma*patch().deltaCoeffs()/Ep);
153 
154  // Restore tag
155  UPstream::msgType() = oldTag;
156 
157  mixedFvPatchScalarField::updateCoeffs();
158 }
159 
160 
162 (
163  Ostream& os
164 ) const
165 {
168  writeEntryIfDifferent<word>(os, "T", "T", TName_);
169 }
170 
171 
172 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
173 
174 namespace Foam
175 {
177  (
180  );
181 }
182 
183 
184 // ************************************************************************* //
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.
A 'mixed' boundary condition that implements a Marshak condition for the incident radiation field (us...
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.
MarshakRadiationFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
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
const Type & value() const
Return const reference to value.
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
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:249
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
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.
A class for handling words, derived from string.
Definition: word.H:62
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m^2/K^4].
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
Namespace for OpenFOAM.
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
dimensionedScalar pow4(const dimensionedScalar &ds)
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
dictionary dict
volScalarField & p