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-2021 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 "fvPatchFieldMapper.H"
29 #include "volFields.H"
30 #include "radiationModel.H"
32 
33 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
34 
37 (
38  const fvPatch& p,
40 )
41 :
42  mixedFvPatchScalarField(p, iF),
43  radiationCoupledBase(p, "undefined", scalarField::null()),
44  Trad_(p.size())
45 {
46  refValue() = 0.0;
47  refGrad() = 0.0;
48  valueFraction() = 0.0;
49 }
50 
51 
54 (
55  const fvPatch& p,
57  const dictionary& dict
58 )
59 :
60  mixedFvPatchScalarField(p, iF),
61  radiationCoupledBase(p, dict),
62  Trad_("Trad", dict, p.size())
63 {
64  // refValue updated on each call to updateCoeffs()
65  refValue() = 4.0*constant::physicoChemical::sigma.value()*pow4(Trad_);
66 
67  // zero gradient
68  refGrad() = 0.0;
69 
70  valueFraction() = 1.0;
71 
72  fvPatchScalarField::operator=(refValue());
73 }
74 
75 
78 (
80  const fvPatch& p,
82  const fvPatchFieldMapper& mapper
83 )
84 :
85  mixedFvPatchScalarField(ptf, p, iF, mapper),
87  (
88  p,
89  ptf.emissivityMethod(),
90  ptf.emissivity_,
91  mapper
92  ),
93  Trad_(mapper(ptf.Trad_))
94 {}
95 
96 
99 (
102 )
103 :
104  mixedFvPatchScalarField(ptf, iF),
106  (
107  ptf.patch(),
108  ptf.emissivityMethod(),
109  ptf.emissivity_
110  ),
111  Trad_(ptf.Trad_)
112 {}
113 
114 
115 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
116 
118 (
119  const fvPatchFieldMapper& m
120 )
121 {
122  mixedFvPatchScalarField::autoMap(m);
123  radiationCoupledBase::autoMap(m);
124  m(Trad_, Trad_);
125 }
126 
127 
129 (
130  const fvPatchScalarField& ptf,
131  const labelList& addr
132 )
133 {
134  mixedFvPatchScalarField::rmap(ptf, addr);
135  radiationCoupledBase::rmap(ptf, addr);
137  refCast<const MarshakRadiationFixedTemperatureFvPatchScalarField>(ptf);
138 
139  Trad_.rmap(mrptf.Trad_, addr);
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  // Re-calc reference value
156  refValue() = 4.0*constant::physicoChemical::sigma.value()*pow4(Trad_);
157 
158  // Diffusion coefficient - created by radiation model's ::updateCoeffs()
159  const scalarField& gamma =
160  patch().lookupPatchField<volScalarField, scalar>("gammaRad");
161 
162  const scalarField temissivity = emissivity();
163 
164  const scalarField Ep(temissivity/(2*(2 - temissivity)));
165 
166  // Set value fraction
167  valueFraction() = 1.0/(1.0 + gamma*patch().deltaCoeffs()/Ep);
168 
169  // Restore tag
170  UPstream::msgType() = oldTag;
171 
172  mixedFvPatchScalarField::updateCoeffs();
173 }
174 
175 
177 (
178  Ostream& os
179 ) const
180 {
183  writeEntry(os, "Trad", Trad_);
184 }
185 
186 
187 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
188 
189 namespace Foam
190 {
192  (
195  );
196 }
197 
198 
199 // ************************************************************************* //
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:62
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
scalarField emissivity() const
Calculate corresponding emissivity field.
static int & msgType()
Message tag of standard messages.
Definition: UPstream.H:476
scalarField emissivity_
Emissivity.
Macros for easy insertion into run-time selection tables.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:57
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m^2/K^4].
MarshakRadiationFixedTemperatureFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
void write(Ostream &) const
Write.
Foam::fvPatchFieldMapper.
const Type & value() const
Return const reference to value.
void write(std::ostream &os, const bool binary, List< floatScalar > &fField)
Write floats ascii or binary.
A &#39;mixed&#39; boundary condition that implements a Marshak condition for the incident radiation field (us...
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Common functions to emissivity. It gets supplied from lookup into a dictionary or calculated by the s...
word emissivityMethod() const
Method to obtain emissivity.
virtual label size() const
Return size.
Definition: fvPatch.H:156
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensionedScalar pow4(const dimensionedScalar &ds)
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Namespace for OpenFOAM.