MarshakRadiationFixedTemperatureFvPatchScalarField.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 
37 (
38  const fvPatch& p,
40  const dictionary& dict
41 )
42 :
43  mixedFvPatchScalarField(p, iF, dict, false),
45  Trad_("Trad", dimTemperature, dict, p.size())
46 {
47  // refValue updated on each call to updateCoeffs()
48  refValue() = 4.0*constant::physicoChemical::sigma.value()*pow4(Trad_);
49 
50  // zero gradient
51  refGrad() = 0.0;
52 
53  valueFraction() = 1.0;
54 
56 }
57 
58 
61 (
63  const fvPatch& p,
65  const fieldMapper& mapper
66 )
67 :
68  mixedFvPatchScalarField(ptf, p, iF, mapper),
70  (
71  p,
72  ptf.emissivityMethod(),
73  ptf.emissivity_,
74  mapper
75  ),
76  Trad_(mapper(ptf.Trad_))
77 {}
78 
79 
82 (
85 )
86 :
87  mixedFvPatchScalarField(ptf, iF),
89  (
90  ptf.patch(),
91  ptf.emissivityMethod(),
92  ptf.emissivity_
93  ),
94  Trad_(ptf.Trad_)
95 {}
96 
97 
98 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
99 
101 (
102  const fvPatchScalarField& ptf,
103  const fieldMapper& mapper
104 )
105 {
106  mixedFvPatchScalarField::map(ptf, mapper);
107  radiationCoupledBase::map(ptf, mapper);
109  refCast<const MarshakRadiationFixedTemperatureFvPatchScalarField>(ptf);
110 
111  mapper(Trad_, mrptf.Trad_);
112 }
113 
114 
116 (
117  const fvPatchScalarField& ptf
118 )
119 {
120  mixedFvPatchScalarField::reset(ptf);
123  refCast<const MarshakRadiationFixedTemperatureFvPatchScalarField>(ptf);
124 
125  Trad_.reset(mrptf.Trad_);
126 }
127 
128 
130 {
131  if (this->updated())
132  {
133  return;
134  }
135 
136  // Since we're inside initEvaluate/evaluate there might be processor
137  // comms underway. Change the tag we use.
138  int oldTag = UPstream::msgType();
139  UPstream::msgType() = oldTag+1;
140 
141  // Re-calc reference value
142  refValue() = 4.0*constant::physicoChemical::sigma.value()*pow4(Trad_);
143 
144  // Diffusion coefficient - created by radiation model's ::updateCoeffs()
145  const scalarField& gamma =
146  patch().lookupPatchField<volScalarField, scalar>("gammaRad");
147 
148  const scalarField emissivity(this->emissivity());
149 
150  const scalarField Ep(emissivity/(2*(2 - emissivity)));
151 
152  // Set value fraction
153  valueFraction() = 1.0/(1.0 + gamma*patch().deltaCoeffs()/Ep);
154 
155  // Restore tag
156  UPstream::msgType() = oldTag;
157 
158  mixedFvPatchScalarField::updateCoeffs();
159 }
160 
161 
163 (
164  Ostream& os
165 ) const
166 {
169  writeEntry(os, "Trad", Trad_);
170 }
171 
172 
173 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
174 
175 namespace Foam
176 {
178  (
181  );
182 }
183 
184 
185 // ************************************************************************* //
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...
A 'mixed' boundary condition that implements a Marshak condition for the incident radiation field (us...
MarshakRadiationFixedTemperatureFvPatchScalarField(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.
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.
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.
const dimensionSet dimTemperature
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
dimensionedScalar pow4(const dimensionedScalar &ds)
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
dictionary dict
volScalarField & p