epsilonWallFunctionFvPatchScalarField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2019 OpenFOAM Foundation
9  Copyright (C) 2017-2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
31 #include "turbulenceModel.H"
32 #include "fvMatrix.H"
34 
35 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
36 
37 const Foam::Enum
38 <
39  Foam::epsilonWallFunctionFvPatchScalarField::blendingType
40 >
41 Foam::epsilonWallFunctionFvPatchScalarField::blendingTypeNames
42 ({
43  { blendingType::STEPWISE , "stepwise" },
44  { blendingType::MAX , "max" },
45  { blendingType::BINOMIAL , "binomial" },
46  { blendingType::EXPONENTIAL, "exponential" }
47 });
48 
50 
51 
52 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
53 
55 {
56  if (master_ != -1)
57  {
58  return;
59  }
60 
61  const volScalarField& epsilon =
62  static_cast<const volScalarField&>(this->internalField());
63 
64  const volScalarField::Boundary& bf = epsilon.boundaryField();
65 
66  label master = -1;
67  forAll(bf, patchi)
68  {
69  if (isA<epsilonWallFunctionFvPatchScalarField>(bf[patchi]))
70  {
72 
73  if (master == -1)
74  {
75  master = patchi;
76  }
77 
78  epf.master() = master;
79  }
80  }
81 }
82 
83 
85 {
86  const volScalarField& epsilon =
87  static_cast<const volScalarField&>(this->internalField());
88 
89  const volScalarField::Boundary& bf = epsilon.boundaryField();
90 
91  const fvMesh& mesh = epsilon.mesh();
92 
93  if (initialised_ && !mesh.changing())
94  {
95  return;
96  }
97 
98  volScalarField weights
99  (
100  IOobject
101  (
102  "weights",
103  mesh.time().timeName(),
104  mesh,
107  false // do not register
108  ),
109  mesh,
111  );
112 
113  DynamicList<label> epsilonPatches(bf.size());
114  forAll(bf, patchi)
115  {
116  if (isA<epsilonWallFunctionFvPatchScalarField>(bf[patchi]))
117  {
118  epsilonPatches.append(patchi);
119 
120  const labelUList& faceCells = bf[patchi].patch().faceCells();
121  for (const auto& faceCell : faceCells)
122  {
123  ++weights[faceCell];
124  }
125  }
126  }
127 
128  cornerWeights_.setSize(bf.size());
129 
130  for (const auto& patchi : epsilonPatches)
131  {
132  const fvPatchScalarField& wf = weights.boundaryField()[patchi];
133  cornerWeights_[patchi] = 1.0/wf.patchInternalField();
134  }
135 
136  G_.setSize(internalField().size(), Zero);
137  epsilon_.setSize(internalField().size(), Zero);
138 
139  initialised_ = true;
140 }
141 
142 
145 (
146  const label patchi
147 )
148 {
149  const volScalarField& epsilon =
150  static_cast<const volScalarField&>(this->internalField());
151 
152  const volScalarField::Boundary& bf = epsilon.boundaryField();
153 
155  refCast<const epsilonWallFunctionFvPatchScalarField>(bf[patchi]);
156 
157  return const_cast<epsilonWallFunctionFvPatchScalarField&>(epf);
158 }
159 
160 
162 (
164  scalarField& G0,
166 )
167 {
168  // Accumulate all of the G and epsilon contributions
169  forAll(cornerWeights_, patchi)
170  {
171  if (!cornerWeights_[patchi].empty())
172  {
173  epsilonWallFunctionFvPatchScalarField& epf = epsilonPatch(patchi);
174 
175  const List<scalar>& w = cornerWeights_[patchi];
176 
177  epf.calculate(turbulence, w, epf.patch(), G0, epsilon0);
178  }
179  }
180 
181  // Apply zero-gradient condition for epsilon
182  forAll(cornerWeights_, patchi)
183  {
184  if (!cornerWeights_[patchi].empty())
185  {
186  epsilonWallFunctionFvPatchScalarField& epf = epsilonPatch(patchi);
187 
188  epf == scalarField(epsilon0, epf.patch().faceCells());
189  }
190  }
191 }
192 
193 
195 (
196  const turbulenceModel& turbModel,
197  const List<scalar>& cornerWeights,
198  const fvPatch& patch,
199  scalarField& G0,
201 )
202 {
203  const label patchi = patch.index();
204 
206  nutWallFunctionFvPatchScalarField::nutw(turbModel, patchi);
207 
208  const scalarField& y = turbModel.y()[patchi];
209 
210  const tmp<scalarField> tnuw = turbModel.nu(patchi);
211  const scalarField& nuw = tnuw();
212 
213  const tmp<volScalarField> tk = turbModel.k();
214  const volScalarField& k = tk();
215 
216  const fvPatchVectorField& Uw = turbModel.U().boundaryField()[patchi];
217 
218  const scalarField magGradUw(mag(Uw.snGrad()));
219 
220  const scalar Cmu25 = pow025(nutw.Cmu());
221  const scalar Cmu75 = pow(nutw.Cmu(), 0.75);
222 
223  // Set epsilon and G
224  forAll(nutw, facei)
225  {
226  const label celli = patch.faceCells()[facei];
227 
228  const scalar yPlus = Cmu25*y[facei]*sqrt(k[celli])/nuw[facei];
229 
230  const scalar w = cornerWeights[facei];
231 
232  scalar epsilonBlended = 0.0;
233 
234  // Contribution from the viscous sublayer
235  const scalar epsilonVis = w*2.0*k[celli]*nuw[facei]/sqr(y[facei]);
236 
237  // Contribution from the inertial sublayer
238  const scalar epsilonLog =
239  w*Cmu75*pow(k[celli], 1.5)/(nutw.kappa()*y[facei]);
240 
241  switch (blending_)
242  {
243  case blendingType::STEPWISE:
244  {
245  if (lowReCorrection_ && yPlus < nutw.yPlusLam())
246  {
247  epsilonBlended = epsilonVis;
248  }
249  else
250  {
251  epsilonBlended = epsilonLog;
252  }
253  break;
254  }
255 
256  case blendingType::MAX:
257  {
258  // (PH:Eq. 27)
259  epsilonBlended = max(epsilonVis, epsilonLog);
260  break;
261  }
262 
263  case blendingType::BINOMIAL:
264  {
265  // (ME:Eqs. 15-16)
266  epsilonBlended =
267  pow
268  (
269  pow(epsilonVis, n_) + pow(epsilonLog, n_),
270  1.0/n_
271  );
272  break;
273  }
274 
275  case blendingType::EXPONENTIAL:
276  {
277  // (PH:p. 193)
278  const scalar Gamma = 0.001*pow4(yPlus)/(1.0 + yPlus);
279  const scalar invGamma = 1.0/(Gamma + ROOTVSMALL);
280  epsilonBlended =
281  epsilonVis*exp(-Gamma) + epsilonLog*exp(-invGamma);
282  break;
283  }
284  }
285 
286  epsilon0[celli] += epsilonBlended;
287 
288  if (!(lowReCorrection_ && yPlus < nutw.yPlusLam()))
289  {
290  G0[celli] +=
291  w
292  *(nutw[facei] + nuw[facei])
293  *magGradUw[facei]
294  *Cmu25*sqrt(k[celli])
295  /(nutw.kappa()*y[facei]);
296  }
297  }
298 }
299 
300 
301 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
302 
305 (
306  const fvPatch& p,
308 )
309 :
311  blending_(blendingType::STEPWISE),
312  n_(2.0),
313  lowReCorrection_(false),
314  initialised_(false),
315  master_(-1),
316  G_(),
317  epsilon_(),
318  cornerWeights_()
319 {}
320 
321 
324 (
326  const fvPatch& p,
328  const fvPatchFieldMapper& mapper
329 )
330 :
331  fixedValueFvPatchField<scalar>(ptf, p, iF, mapper),
332  blending_(ptf.blending_),
333  n_(ptf.n_),
334  lowReCorrection_(ptf.lowReCorrection_),
335  initialised_(false),
336  master_(-1),
337  G_(),
338  epsilon_(),
339  cornerWeights_()
340 {}
341 
342 
345 (
346  const fvPatch& p,
348  const dictionary& dict
349 )
350 :
352  blending_
353  (
354  blendingTypeNames.getOrDefault
355  (
356  "blending",
357  dict,
358  blendingType::STEPWISE
359  )
360  ),
361  n_
362  (
363  dict.getCheckOrDefault<scalar>
364  (
365  "n",
366  2.0,
368  )
369  ),
370  lowReCorrection_(dict.getOrDefault("lowReCorrection", false)),
371  initialised_(false),
372  master_(-1),
373  G_(),
374  epsilon_(),
375  cornerWeights_()
376 {
377  // Apply zero-gradient condition on start-up
378  this->operator==(patchInternalField());
379 }
380 
381 
384 (
386 )
387 :
389  blending_(ewfpsf.blending_),
390  n_(ewfpsf.n_),
391  lowReCorrection_(ewfpsf.lowReCorrection_),
392  initialised_(false),
393  master_(-1),
394  G_(),
395  epsilon_(),
396  cornerWeights_()
397 {}
398 
399 
402 (
405 )
406 :
407  fixedValueFvPatchField<scalar>(ewfpsf, iF),
408  blending_(ewfpsf.blending_),
409  n_(ewfpsf.n_),
410  lowReCorrection_(ewfpsf.lowReCorrection_),
411  initialised_(false),
412  master_(-1),
413  G_(),
414  epsilon_(),
415  cornerWeights_()
416 {}
417 
418 
419 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
420 
422 (
423  bool init
424 )
425 {
426  if (patch().index() == master_)
427  {
428  if (init)
429  {
430  G_ = 0.0;
431  }
432 
433  return G_;
434  }
435 
436  return epsilonPatch(master_).G();
437 }
438 
439 
441 (
442  bool init
443 )
444 {
445  if (patch().index() == master_)
446  {
447  if (init)
448  {
449  epsilon_ = 0.0;
450  }
451 
452  return epsilon_;
453  }
454 
455  return epsilonPatch(master_).epsilon(init);
456 }
457 
458 
460 {
461  if (updated())
462  {
463  return;
464  }
465 
466  const turbulenceModel& turbModel = db().lookupObject<turbulenceModel>
467  (
469  (
471  internalField().group()
472  )
473  );
474 
475  setMaster();
476 
477  if (patch().index() == master_)
478  {
479  createAveragingWeights();
480  calculateTurbulenceFields(turbModel, G(true), epsilon(true));
481  }
482 
483  const scalarField& G0 = this->G();
484  const scalarField& epsilon0 = this->epsilon();
485 
486  typedef DimensionedField<scalar, volMesh> FieldType;
487 
488  FieldType& G = db().lookupObjectRef<FieldType>(turbModel.GName());
489 
490  FieldType& epsilon = const_cast<FieldType&>(internalField());
491 
492  forAll(*this, facei)
493  {
494  const label celli = patch().faceCells()[facei];
495 
496  G[celli] = G0[celli];
497  epsilon[celli] = epsilon0[celli];
498  }
499 
501 }
502 
503 
505 (
506  const scalarField& weights
507 )
508 {
509  if (updated())
510  {
511  return;
512  }
513 
514  const turbulenceModel& turbModel = db().lookupObject<turbulenceModel>
515  (
517  (
519  internalField().group()
520  )
521  );
522 
523  setMaster();
524 
525  if (patch().index() == master_)
526  {
527  createAveragingWeights();
528  calculateTurbulenceFields(turbModel, G(true), epsilon(true));
529  }
530 
531  const scalarField& G0 = this->G();
532  const scalarField& epsilon0 = this->epsilon();
533 
534  typedef DimensionedField<scalar, volMesh> FieldType;
535 
536  FieldType& G = db().lookupObjectRef<FieldType>(turbModel.GName());
537 
538  FieldType& epsilon = const_cast<FieldType&>(internalField());
539 
540  scalarField& epsilonf = *this;
541 
542  // Only set the values if the weights are > tolerance
543  forAll(weights, facei)
544  {
545  const scalar w = weights[facei];
546 
547  if (w > tolerance_)
548  {
549  const label celli = patch().faceCells()[facei];
550 
551  G[celli] = (1.0 - w)*G[celli] + w*G0[celli];
552  epsilon[celli] = (1.0 - w)*epsilon[celli] + w*epsilon0[celli];
553  epsilonf[facei] = epsilon[celli];
554  }
555  }
556 
558 }
559 
560 
562 (
563  fvMatrix<scalar>& matrix
564 )
565 {
566  if (manipulatedMatrix())
567  {
568  return;
569  }
570 
571  matrix.setValues(patch().faceCells(), patchInternalField());
572 
574 }
575 
576 
578 (
579  fvMatrix<scalar>& matrix,
580  const Field<scalar>& weights
581 )
582 {
583  if (manipulatedMatrix())
584  {
585  return;
586  }
587 
588  DynamicList<label> constraintCells(weights.size());
589  DynamicList<scalar> constraintValues(weights.size());
590  const labelUList& faceCells = patch().faceCells();
591 
592  const DimensionedField<scalar, volMesh>& fld = internalField();
593 
594  forAll(weights, facei)
595  {
596  // Only set the values if the weights are > tolerance
597  if (weights[facei] > tolerance_)
598  {
599  const label celli = faceCells[facei];
600 
601  constraintCells.append(celli);
602  constraintValues.append(fld[celli]);
603  }
604  }
605 
606  if (debug)
607  {
608  Pout<< "Patch: " << patch().name()
609  << ": number of constrained cells = " << constraintCells.size()
610  << " out of " << patch().size()
611  << endl;
612  }
613 
614  matrix.setValues(constraintCells, constraintValues);
615 
617 }
618 
619 
621 (
622  Ostream& os
623 ) const
624 {
625  os.writeEntry("lowReCorrection", lowReCorrection_);
626  os.writeEntry("blending", blendingTypeNames[blending_]);
627  os.writeEntry("n", n_);
629 }
630 
631 
632 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
633 
634 namespace Foam
635 {
637  (
640  );
641 }
642 
643 
644 // ************************************************************************* //
Foam::expressions::patchExpr::debug
int debug
Static debugging option.
Foam::fvPatchField< scalar >
Foam::IOobject::NO_WRITE
Definition: IOobject.H:195
Foam::fvPatchField< scalar >::internalField
const DimensionedField< scalar, volMesh > & internalField() const
Return dimensioned internal field reference.
Definition: fvPatchField.H:363
Foam::fvPatchField::snGrad
virtual tmp< Field< Type > > snGrad() const
Return patch-normal gradient.
Definition: fvPatchField.C:225
Foam::scalarField
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
Definition: primitiveFieldsFwd.H:52
Foam::IOobject
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:169
Foam::turbulenceModel::nu
virtual tmp< volScalarField > nu() const =0
Return the laminar viscosity.
p
volScalarField & p
Definition: createFieldRefs.H:8
Foam::Enum
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: IOstreamOption.H:57
Foam::epsilonWallFunctionFvPatchScalarField::master
virtual label & master()
Return non-const access to the master patch ID.
Definition: epsilonWallFunctionFvPatchScalarField.H:349
Foam::nutWallFunctionFvPatchScalarField::Cmu
scalar Cmu() const
Return Cmu.
Definition: nutWallFunctionFvPatchScalarField.H:369
turbulence
Info<< "Reading field U\n"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);volScalarField rho(IOobject("rho", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::AUTO_WRITE), thermo.rho());volVectorField rhoU(IOobject("rhoU", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *U);volScalarField rhoE(IOobject("rhoE", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE), rho *(e+0.5 *magSqr(U)));surfaceScalarField pos(IOobject("pos", runTime.timeName(), mesh), mesh, dimensionedScalar("pos", dimless, 1.0));surfaceScalarField neg(IOobject("neg", runTime.timeName(), mesh), mesh, dimensionedScalar("neg", dimless, -1.0));surfaceScalarField phi("phi", fvc::flux(rhoU));Info<< "Creating turbulence model\n"<< endl;autoPtr< compressible::turbulenceModel > turbulence(compressible::turbulenceModel::New(rho, U, phi, thermo))
Definition: createFields.H:94
Foam::constant::atomic::group
constexpr const char *const group
Group name for atomic constants.
Definition: atomicConstants.H:52
Foam::constant::universal::G
const dimensionedScalar G
Newtonian constant of gravitation.
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::DynamicList< label >
Foam::epsilonWallFunctionFvPatchScalarField::updateWeightedCoeffs
virtual void updateWeightedCoeffs(const scalarField &weights)
Update the coefficients associated with the patch field.
Definition: epsilonWallFunctionFvPatchScalarField.C:505
Foam::epsilonWallFunctionFvPatchScalarField::tolerance_
static scalar tolerance_
Tolerance used in weighted calculations.
Definition: epsilonWallFunctionFvPatchScalarField.H:293
Foam::Time::timeName
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
Foam::turbulenceModel::propertiesName
static const word propertiesName
Default name of the turbulence properties dictionary.
Definition: turbulenceModel.H:100
nutWallFunctionFvPatchScalarField.H
fvMatrix.H
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::exp
dimensionedScalar exp(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:261
epsilonWallFunctionFvPatchScalarField.H
Foam::Pout
prefixOSstream Pout
OSstream wrapped stdout (std::cout) with parallel prefix.
Foam::turbulenceModel::y
const nearWallDist & y() const
Return the near wall distances.
Definition: turbulenceModel.H:159
Foam::epsilonWallFunctionFvPatchScalarField::createAveragingWeights
virtual void createAveragingWeights()
Definition: epsilonWallFunctionFvPatchScalarField.C:84
Foam::epsilonWallFunctionFvPatchScalarField::epsilonWallFunctionFvPatchScalarField
epsilonWallFunctionFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &)
Construct from patch and internal field.
Definition: epsilonWallFunctionFvPatchScalarField.C:305
Foam::fixedValueFvPatchField< scalar >
Foam::constant::electromagnetic::epsilon0
const dimensionedScalar epsilon0
Electric constant: default SI units: [F/m].
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::fixedValueFvPatchField::write
virtual void write(Ostream &) const
Write.
Definition: fixedValueFvPatchField.C:156
Foam::pow025
dimensionedScalar pow025(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:133
Foam::pow4
dimensionedScalar pow4(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:100
Foam::epsilonWallFunctionFvPatchScalarField::epsilon
scalarField & epsilon(bool init=false)
Return non-const access to the master's epsilon field.
Definition: epsilonWallFunctionFvPatchScalarField.C:441
Foam::Field< scalar >
Foam::operator==
tmp< faMatrix< Type > > operator==(const faMatrix< Type > &, const faMatrix< Type > &)
Foam::constant::electromagnetic::G0
const dimensionedScalar G0
Conductance quantum: default SI units: [S].
Foam::fvPatch
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:65
Foam::epsilonWallFunctionFvPatchScalarField
This boundary condition provides a wall constraint on the turbulent kinetic energy dissipation rate,...
Definition: epsilonWallFunctionFvPatchScalarField.H:259
init
mesh init(true)
Foam::fvPatchField::patchInternalField
virtual tmp< Field< Type > > patchInternalField() const
Return internal field next to patch as patch field.
Definition: fvPatchField.C:233
Foam::epsilonWallFunctionFvPatchScalarField::calculate
virtual void calculate(const turbulenceModel &turbulence, const List< scalar > &cornerWeights, const fvPatch &patch, scalarField &G, scalarField &epsilon)
Calculate the epsilon and G.
Definition: epsilonWallFunctionFvPatchScalarField.C:195
Foam::epsilonWallFunctionFvPatchScalarField::lowReCorrection_
const bool lowReCorrection_
Apply low-Re correction term (default = no)
Definition: epsilonWallFunctionFvPatchScalarField.H:296
Foam::dimensionedScalar
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
Definition: dimensionedScalarFwd.H:42
Foam::epsilonWallFunctionFvPatchScalarField::write
virtual void write(Ostream &) const
Write.
Definition: epsilonWallFunctionFvPatchScalarField.C:621
Foam::dictionary::getCheckOrDefault
T getCheckOrDefault(const word &keyword, const T &deflt, const Predicate &pred, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:209
Foam::pow
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Definition: dimensionedScalar.C:75
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
fld
gmvFile<< "tracers "<< particles.size()<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().x()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().y()<< ' ';}gmvFile<< nl;for(const passiveParticle &p :particles){ gmvFile<< p.position().z()<< ' ';}gmvFile<< nl;for(const word &name :lagrangianScalarNames){ IOField< scalar > fld(IOobject(name, runTime.timeName(), cloud::prefix, mesh, IOobject::MUST_READ, IOobject::NO_WRITE))
Definition: gmvOutputLagrangian.H:23
Foam::turbulenceModel::GName
word GName() const
Helper function to return the name of the turbulence G field.
Definition: turbulenceModel.H:138
Foam::epsilonWallFunctionFvPatchScalarField::G
scalarField & G(bool init=false)
Return non-const access to the master's G field.
Definition: epsilonWallFunctionFvPatchScalarField.C:422
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::turbulenceModel
Abstract base class for turbulence models (RAS, LES and laminar).
Definition: turbulenceModel.H:63
Foam::epsilonWallFunctionFvPatchScalarField::updateCoeffs
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: epsilonWallFunctionFvPatchScalarField.C:459
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
os
OBJstream os(runTime.globalPath()/outputName)
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam::turbulenceModel::k
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
Foam::fvMesh
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:85
Foam::epsilonWallFunctionFvPatchScalarField::manipulateMatrix
virtual void manipulateMatrix(fvMatrix< scalar > &matrix)
Manipulate matrix.
Definition: epsilonWallFunctionFvPatchScalarField.C:562
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::nutWallFunctionFvPatchScalarField::yPlusLam
static scalar yPlusLam(const scalar kappa, const scalar E)
Estimate the y+ at the intersection of the two sublayers.
Definition: nutWallFunctionFvPatchScalarField.C:250
Foam::fvPatchField::manipulateMatrix
virtual void manipulateMatrix(fvMatrix< Type > &matrix)
Manipulate matrix.
Definition: fvPatchField.C:354
Foam::epsilonWallFunctionFvPatchScalarField::calculateTurbulenceFields
virtual void calculateTurbulenceFields(const turbulenceModel &turbulence, scalarField &G0, scalarField &epsilon0)
Main driver to calculate the turbulence fields.
Definition: epsilonWallFunctionFvPatchScalarField.C:162
Foam::fvPatchField::updateCoeffs
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Definition: fvPatchField.C:321
Foam::epsilonWallFunctionFvPatchScalarField::epsilonPatch
virtual epsilonWallFunctionFvPatchScalarField & epsilonPatch(const label patchi)
Helper function to return non-const access to an epsilon patch.
Definition: epsilonWallFunctionFvPatchScalarField.C:145
Foam::fvMatrix::setValues
void setValues(const labelUList &cellLabels, const Type &value)
Definition: fvMatrix.C:1059
Foam::polyMesh::changing
bool changing() const noexcept
Is mesh changing (topology changing and/or moving)
Definition: polyMesh.H:550
Foam::sqr
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Definition: dimensionedSymmTensor.C:51
Foam::foamVersion::patch
const std::string patch
OpenFOAM patch number as a std::string.
Foam::MinMax::ge
static MinMax< T > ge(const T &minVal)
A semi-infinite range from minVal to the type max.
Definition: MinMaxI.H:31
Foam::nutWallFunctionFvPatchScalarField::nutw
static const nutWallFunctionFvPatchScalarField & nutw(const turbulenceModel &turbModel, const label patchi)
Return the nut patchField for the given wall patch.
Definition: nutWallFunctionFvPatchScalarField.C:235
Foam::List< scalar >
Foam::sqrt
dimensionedScalar sqrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:144
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::nutWallFunctionFvPatchScalarField::kappa
scalar kappa() const
Return kappa.
Definition: nutWallFunctionFvPatchScalarField.H:375
Foam::UList< label >
k
label k
Boltzmann constant.
Definition: LISASMDCalcMethod2.H:41
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
Foam::epsilonWallFunctionFvPatchScalarField::master_
label master_
Master patch ID.
Definition: epsilonWallFunctionFvPatchScalarField.H:302
Foam::fvMatrix
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvPatchField.H:68
Foam::Ostream::writeEntry
Ostream & writeEntry(const keyType &key, const T &value)
Write a keyword/value entry.
Definition: Ostream.H:236
Foam::epsilonWallFunctionFvPatchScalarField::setMaster
virtual void setMaster()
Definition: epsilonWallFunctionFvPatchScalarField.C:54
Foam::fvPatch::faceCells
virtual const labelUList & faceCells() const
Return faceCells.
Definition: fvPatch.C:113
Foam::fvMesh::time
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:280
yPlus
scalar yPlus
Definition: evaluateNearWall.H:16
Foam::fvPatchFieldMapper
Foam::fvPatchFieldMapper.
Definition: fvPatchFieldMapper.H:47
Foam::fvPatchField::patch
const fvPatch & patch() const
Return patch.
Definition: fvPatchField.H:357
epsilon
scalar epsilon
Definition: evaluateNearWall.H:7
Foam::IOobject::groupName
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::makePatchTypeField
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Foam::dictionary::getOrDefault
T getOrDefault(const word &keyword, const T &deflt, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:148
Foam::GeometricField< scalar, fvPatchField, volMesh >
Foam::nutWallFunctionFvPatchScalarField
The class nutWallFunction is a base class that parents the derived boundary conditions which provide ...
Definition: nutWallFunctionFvPatchScalarField.H:251
Foam::IOobject::NO_READ
Definition: IOobject.H:188
Foam::faceCells
Smooth ATC in cells next to a set of patches supplied by type.
Definition: faceCells.H:56
turbulenceModel.H
Foam::GeometricField::boundaryField
const Boundary & boundaryField() const
Return const-reference to the boundary field.
Definition: GeometricFieldI.H:62
Foam::DimensionedField
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
Definition: DimensionedField.H:54
Foam::dimless
const dimensionSet dimless
Dimensionless.
Definition: dimensionSets.C:189
Foam::turbulenceModel::U
const volVectorField & U() const
Access function to velocity field.
Definition: turbulenceModel.H:144
y
scalar y
Definition: LISASMDCalcMethod1.H:14