alphatJayatillekeWallFunctionFvPatchScalarField.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-2018 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"
31 #include "wallFvPatch.H"
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 namespace compressible
38 {
39 
40 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
41 
42 scalar alphatJayatillekeWallFunctionFvPatchScalarField::maxExp_ = 50.0;
43 scalar alphatJayatillekeWallFunctionFvPatchScalarField::tolerance_ = 0.01;
44 label alphatJayatillekeWallFunctionFvPatchScalarField::maxIters_ = 10;
45 
46 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
47 
48 void alphatJayatillekeWallFunctionFvPatchScalarField::checkType()
49 {
50  if (!isA<wallFvPatch>(patch()))
51  {
53  << "Patch type for patch " << patch().name() << " must be wall\n"
54  << "Current patch type is " << patch().type() << nl
55  << exit(FatalError);
56  }
57 }
58 
59 
60 scalar alphatJayatillekeWallFunctionFvPatchScalarField::Psmooth
61 (
62  const scalar Prat
63 ) const
64 {
65  return 9.24*(pow(Prat, 0.75) - 1.0)*(1.0 + 0.28*exp(-0.007*Prat));
66 }
67 
68 
69 scalar alphatJayatillekeWallFunctionFvPatchScalarField::yPlusTherm
70 (
71  const scalar P,
72  const scalar Prat
73 ) const
74 {
75  scalar ypt = 11.0;
76 
77  for (int i=0; i<maxIters_; i++)
78  {
79  scalar f = ypt - (log(E_*ypt)/kappa_ + P)/Prat;
80  scalar df = 1.0 - 1.0/(ypt*kappa_*Prat);
81  scalar yptNew = ypt - f/df;
82 
83  if (yptNew < vSmall)
84  {
85  return 0;
86  }
87  else if (mag(yptNew - ypt) < tolerance_)
88  {
89  return yptNew;
90  }
91  else
92  {
93  ypt = yptNew;
94  }
95  }
96 
97  return ypt;
98 }
99 
100 
101 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
102 
105 (
106  const fvPatch& p,
108 )
109 :
110  fixedValueFvPatchScalarField(p, iF),
111  Prt_(0.85),
112  Cmu_(0.09),
113  kappa_(0.41),
114  E_(9.8)
115 {
116  checkType();
117 }
118 
119 
122 (
124  const fvPatch& p,
126  const fvPatchFieldMapper& mapper
127 )
128 :
129  fixedValueFvPatchScalarField(ptf, p, iF, mapper),
130  Prt_(ptf.Prt_),
131  Cmu_(ptf.Cmu_),
132  kappa_(ptf.kappa_),
133  E_(ptf.E_)
134 {}
135 
136 
139 (
140  const fvPatch& p,
142  const dictionary& dict
143 )
144 :
145  fixedValueFvPatchScalarField(p, iF, dict),
146  Prt_(dict.lookupOrDefault<scalar>("Prt", 0.85)),
147  Cmu_(dict.lookupOrDefault<scalar>("Cmu", 0.09)),
148  kappa_(dict.lookupOrDefault<scalar>("kappa", 0.41)),
149  E_(dict.lookupOrDefault<scalar>("E", 9.8))
150 {
151  checkType();
152 }
153 
154 
157 (
159 )
160 :
161  fixedValueFvPatchScalarField(awfpsf),
162  Prt_(awfpsf.Prt_),
163  Cmu_(awfpsf.Cmu_),
164  kappa_(awfpsf.kappa_),
165  E_(awfpsf.E_)
166 {
167  checkType();
168 }
169 
170 
173 (
176 )
177 :
178  fixedValueFvPatchScalarField(awfpsf, iF),
179  Prt_(awfpsf.Prt_),
180  Cmu_(awfpsf.Cmu_),
181  kappa_(awfpsf.kappa_),
182  E_(awfpsf.E_)
183 {
184  checkType();
185 }
186 
187 
188 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
189 
191 {
192  if (updated())
193  {
194  return;
195  }
196 
197  const label patchi = patch().index();
198 
199  // Retrieve turbulence properties from model
200  const compressible::turbulenceModel& turbModel =
201  db().lookupObject<compressible::turbulenceModel>
202  (
204  (
205  compressible::turbulenceModel::propertiesName,
206  internalField().group()
207  )
208  );
209 
210  const scalar Cmu25 = pow025(Cmu_);
211 
212  const scalarField& y = turbModel.y()[patchi];
213 
214  const tmp<scalarField> tmuw = turbModel.mu(patchi);
215  const scalarField& muw = tmuw();
216 
217  const tmp<scalarField> talphaw = turbModel.alpha(patchi);
218  const scalarField& alphaw = talphaw();
219 
220  scalarField& alphatw = *this;
221 
222  const tmp<volScalarField> tk = turbModel.k();
223  const volScalarField& k = tk();
224 
225  const fvPatchVectorField& Uw = turbModel.U().boundaryField()[patchi];
226  const scalarField magUp(mag(Uw.patchInternalField() - Uw));
227  const scalarField magGradUw(mag(Uw.snGrad()));
228 
229  const scalarField& rhow = turbModel.rho().boundaryField()[patchi];
230  const fvPatchScalarField& hew =
231  turbModel.transport().he().boundaryField()[patchi];
232 
233  // Heat flux [W/m2] - lagging alphatw
234  const scalarField qDot
235  (
236  turbModel.transport().alphaEff(alphatw, patchi)*hew.snGrad()
237  );
238 
239  // Populate boundary values
240  forAll(alphatw, facei)
241  {
242  label celli = patch().faceCells()[facei];
243 
244  scalar uTau = Cmu25*sqrt(k[celli]);
245 
246  scalar yPlus = uTau*y[facei]/(muw[facei]/rhow[facei]);
247 
248  // Molecular Prandtl number
249  scalar Pr = muw[facei]/alphaw[facei];
250 
251  // Molecular-to-turbulent Prandtl number ratio
252  scalar Prat = Pr/Prt_;
253 
254  // Thermal sublayer thickness
255  scalar P = Psmooth(Prat);
256  scalar yPlusTherm = this->yPlusTherm(P, Prat);
257 
258  // Evaluate new effective thermal diffusivity
259  scalar alphaEff = 0.0;
260  if (yPlus < yPlusTherm)
261  {
262  scalar A = qDot[facei]*rhow[facei]*uTau*y[facei];
263  scalar B = qDot[facei]*Pr*yPlus;
264  scalar C = Pr*0.5*rhow[facei]*uTau*sqr(magUp[facei]);
265  alphaEff = A/(B + C + vSmall);
266  }
267  else
268  {
269  scalar A = qDot[facei]*rhow[facei]*uTau*y[facei];
270  scalar B = qDot[facei]*Prt_*(1.0/kappa_*log(E_*yPlus) + P);
271  scalar magUc = uTau/kappa_*log(E_*yPlusTherm) - mag(Uw[facei]);
272  scalar C =
273  0.5*rhow[facei]*uTau
274  *(Prt_*sqr(magUp[facei]) + (Pr - Prt_)*sqr(magUc));
275  alphaEff = A/(B + C + vSmall);
276  }
277 
278  // Update turbulent thermal diffusivity
279  alphatw[facei] = max(0.0, alphaEff - alphaw[facei]);
280 
281  if (debug)
282  {
283  Info<< " uTau = " << uTau << nl
284  << " Pr = " << Pr << nl
285  << " Prt = " << Prt_ << nl
286  << " qDot = " << qDot[facei] << nl
287  << " yPlus = " << yPlus << nl
288  << " yPlusTherm = " << yPlusTherm << nl
289  << " alphaEff = " << alphaEff << nl
290  << " alphaw = " << alphaw[facei] << nl
291  << " alphatw = " << alphatw[facei] << nl
292  << endl;
293  }
294  }
295 
297 }
298 
299 
301 {
303  os.writeKeyword("Prt") << Prt_ << token::END_STATEMENT << nl;
304  os.writeKeyword("Cmu") << Cmu_ << token::END_STATEMENT << nl;
305  os.writeKeyword("kappa") << kappa_ << token::END_STATEMENT << nl;
306  os.writeKeyword("E") << E_ << token::END_STATEMENT << nl;
307  writeEntry("value", os);
308 }
309 
310 
311 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
312 
314 (
317 );
318 
319 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
320 
321 } // End namespace compressible
322 } // End namespace Foam
323 
324 // ************************************************************************* //
const char *const group
Group name for atomic constants.
dimensionedScalar Pr("Pr", dimless, laminarTransport)
Graphite solid properties.
Definition: C.H:48
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:216
dimensionedScalar log(const dimensionedScalar &ds)
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
dimensionedSymmTensor sqr(const dimensionedVector &dv)
This boundary condition provides a thermal wall function for turbulent thermal diffusivity (usuallyal...
dimensionedScalar sqrt(const dimensionedScalar &ds)
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
dimensionedScalar pow025(const dimensionedScalar &ds)
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
label k
Boltzmann constant.
Templated wrapper class to provide compressible turbulence models thermal diffusivity based thermal t...
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:362
Macros for easy insertion into run-time selection tables.
scalar magUp
virtual tmp< volScalarField > alphaEff() const
Effective thermal turbulent diffusivity of mixture [kg/m/s].
scalar uTau
alphatJayatillekeWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
scalar y
dimensionedScalar exp(const dimensionedScalar &ds)
static word groupName(Name name, const word &group)
Foam::fvPatchFieldMapper.
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:53
static const char nl
Definition: Ostream.H:265
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:224
Ostream & writeKeyword(const keyType &)
Write the keyword followed by an appropriate indentation.
Definition: Ostream.C:54
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
label patchi
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
scalar yPlus
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:311
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
messageStream Info
dimensioned< scalar > mag(const dimensioned< Type > &)
volScalarField alphaEff("alphaEff", turbulence->nu()/Pr+alphat)
A class for managing temporary objects.
Definition: PtrList.H:53
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
virtual tmp< volScalarField > alpha() const
Thermal diffusivity for enthalpy of mixture [kg/m/s].
makePatchTypeField(fvPatchScalarField, thermalBaffleFvPatchScalarField)
Namespace for OpenFOAM.