turbulentMixingLengthFrequencyInletFvPatchScalarField.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 "surfaceFields.H"
30 #include "volFields.H"
31 #include "momentumTransportModel.H"
32 
33 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
34 
35 namespace Foam
36 {
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
42 (
43  const fvPatch& p,
45  const dictionary& dict
46 )
47 :
48  inletOutletFvPatchScalarField(p, iF),
49  mixingLength_(dict.lookup<scalar>("mixingLength", dimLength)),
50  kName_(dict.lookupOrDefault<word>("k", "k"))
51 {
52  this->phiName_ = dict.lookupOrDefault<word>("phi", "phi");
53 
55  (
56  scalarField("value", iF.dimensions(), dict, p.size())
57  );
58 
59  this->refValue() = 0.0;
60  this->refGrad() = 0.0;
61  this->valueFraction() = 0.0;
62 }
63 
64 
67 (
69  const fvPatch& p,
71  const fieldMapper& mapper
72 )
73 :
74  inletOutletFvPatchScalarField(ptf, p, iF, mapper),
75  mixingLength_(ptf.mixingLength_),
76  kName_(ptf.kName_)
77 {}
78 
79 
82 (
85 )
86 :
87  inletOutletFvPatchScalarField(ptf, iF),
88  mixingLength_(ptf.mixingLength_),
89  kName_(ptf.kName_)
90 {}
91 
92 
93 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
94 
96 {
97  if (updated())
98  {
99  return;
100  }
101 
102  // Lookup Cmu corresponding to the turbulence model selected
103 
104  const momentumTransportModel& turbModel =
105  db().lookupType<momentumTransportModel>(internalField().group());
106 
107  const scalar Cmu =
108  turbModel.coeffDict().lookupOrDefault<scalar>("Cmu", 0.09);
109 
110  const scalar Cmu25 = pow(Cmu, 0.25);
111 
112  const fvPatchScalarField& kp =
113  patch().lookupPatchField<volScalarField, scalar>(kName_);
114 
115  const fvsPatchScalarField& phip =
116  patch().lookupPatchField<surfaceScalarField, scalar>(this->phiName_);
117 
118  this->refValue() = sqrt(kp)/(Cmu25*mixingLength_);
119  this->valueFraction() = neg(phip);
120 
121  inletOutletFvPatchScalarField::updateCoeffs();
122 }
123 
124 
126 (
127  Ostream& os
128 ) const
129 {
131  writeEntry(os, "mixingLength", mixingLength_);
132  writeEntry(os, "phi", this->phiName_);
133  writeEntry(os, "k", kName_);
134  writeEntry(os, "value", *this);
135 }
136 
137 
138 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
139 
141 (
144 );
145 
146 
147 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
148 
149 } // End namespace Foam
150 
151 // ************************************************************************* //
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...
const dimensionSet & dimensions() const
Return dimensions.
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T, if not found return the given default.
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 write(Ostream &) const
Write.
Definition: fvPatchField.C:229
friend Ostream & operator(Ostream &, const fvPatchField< Type > &)
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:83
Abstract base class for turbulence models (RAS, LES and laminar).
virtual const dictionary & coeffDict() const =0
Const access to the coefficients dictionary.
This boundary condition provides a turbulence specific dissipation, (omega) inlet condition based on...
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
turbulentMixingLengthFrequencyInletFvPatchScalarField(const fvPatch &, const DimensionedField< scalar, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
A class for handling words, derived from string.
Definition: word.H:62
Namespace for OpenFOAM.
SurfaceField< scalar > surfaceScalarField
const dimensionSet dimLength
Field< scalar > scalarField
Specialisation of Field<T> for scalar.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
dimensionedScalar sqrt(const dimensionedScalar &ds)
VolField< scalar > volScalarField
Definition: volFieldsFwd.H:64
dimensionedScalar neg(const dimensionedScalar &ds)
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
dictionary dict
volScalarField & p
Foam::surfaceFields.