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-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 
36 (
37  const fvPatch& p,
39 )
40 :
41  mixedFvPatchScalarField(p, iF),
42  radiationCoupledBase(p, "undefined", scalarField::null()),
43  TName_("T")
44 {
45  refValue() = 0.0;
46  refGrad() = 0.0;
47  valueFraction() = 0.0;
48 }
49 
50 
52 (
53  const fvPatch& p,
55  const dictionary& dict
56 )
57 :
58  mixedFvPatchScalarField(p, iF),
59  radiationCoupledBase(p, dict),
60  TName_(dict.lookupOrDefault<word>("T", "T"))
61 {
62  if (dict.found("value"))
63  {
64  refValue() = scalarField("value", dict, p.size());
65  }
66  else
67  {
68  refValue() = 0.0;
69  }
70 
71  // zero gradient
72  refGrad() = 0.0;
73 
74  valueFraction() = 1.0;
75 
76  fvPatchScalarField::operator=(refValue());
77 }
78 
79 
81 (
83  const fvPatch& p,
85  const fvPatchFieldMapper& mapper
86 )
87 :
88  mixedFvPatchScalarField(ptf, p, iF, mapper),
90  (
91  p,
92  ptf.emissivityMethod(),
93  ptf.emissivity_,
94  mapper
95  ),
96  TName_(ptf.TName_)
97 {}
98 
99 
101 (
104 )
105 :
106  mixedFvPatchScalarField(ptf, iF),
108  (
109  ptf.patch(),
110  ptf.emissivityMethod(),
111  ptf.emissivity_
112  ),
113  TName_(ptf.TName_)
114 {}
115 
116 
117 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
118 
120 (
121  const fvPatchFieldMapper& m
122 )
123 {
124  mixedFvPatchScalarField::autoMap(m);
125  radiationCoupledBase::autoMap(m);
126 }
127 
128 
130 (
131  const fvPatchScalarField& ptf,
132  const labelList& addr
133 )
134 {
135  mixedFvPatchScalarField::rmap(ptf, addr);
136  radiationCoupledBase::rmap(ptf, addr);
137 }
138 
139 
141 {
142  if (this->updated())
143  {
144  return;
145  }
146 
147  // Since we're inside initEvaluate/evaluate there might be processor
148  // comms underway. Change the tag we use.
149  int oldTag = UPstream::msgType();
150  UPstream::msgType() = oldTag+1;
151 
152  // Temperature field
153  const scalarField& Tp =
154  patch().lookupPatchField<volScalarField, scalar>(TName_);
155 
156  // Re-calc reference value
157  refValue() = 4.0*constant::physicoChemical::sigma.value()*pow4(Tp);
158 
159  // Diffusion coefficient - created by radiation model's ::updateCoeffs()
160  const scalarField& gamma =
161  patch().lookupPatchField<volScalarField, scalar>("gammaRad");
162 
163  const scalarField temissivity = emissivity();
164 
165  const scalarField Ep(temissivity/(2.0*(2.0 - temissivity)));
166 
167  // Set value fraction
168  valueFraction() = 1.0/(1.0 + gamma*patch().deltaCoeffs()/Ep);
169 
170  // Restore tag
171  UPstream::msgType() = oldTag;
172 
173  mixedFvPatchScalarField::updateCoeffs();
174 }
175 
176 
178 (
179  Ostream& os
180 ) const
181 {
184  writeEntryIfDifferent<word>(os, "T", "T", TName_);
185 }
186 
187 
188 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
189 
190 namespace Foam
191 {
193  (
196  );
197 }
198 
199 
200 // ************************************************************************* //
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:643
virtual void autoMap(const fvPatchFieldMapper &)
Map (and resize as needed) from self given a mapping object.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
virtual void rmap(const fvPatchScalarField &, const labelList &)
Reverse map the given fvPatchField onto this fvPatchField.
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].
A class for handling words, derived from string.
Definition: word.H:59
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.
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
volScalarField scalarField(fieldObject, mesh)
MarshakRadiationFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensionedScalar pow4(const dimensionedScalar &ds)
A &#39;mixed&#39; boundary condition that implements a Marshak condition for the incident radiation field (us...
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Namespace for OpenFOAM.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.