multiphaseStabilizedTurbulence.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) 2019-2021 OpenCFD Ltd.
9-------------------------------------------------------------------------------
10License
11 This file is part of OpenFOAM.
12
13 OpenFOAM is free software: you can redistribute it and/or modify it
14 under the terms of the GNU General Public License as published by
15 the Free Software Foundation, either version 3 of the License, or
16 (at your option) any later version.
17
18 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21 for more details.
22
23 You should have received a copy of the GNU General Public License
24 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25
26\*---------------------------------------------------------------------------*/
27
29#include "fvMatrices.H"
31#include "gravityMeshObject.H"
33
34// * * * * * * * * * * * * * Static Member Functions * * * * * * * * * * * * //
35
36namespace Foam
37{
38namespace fv
39{
42 (
43 option,
46 );
47}
48}
49
50// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51
53(
54 const word& sourceName,
55 const word& modelType,
56 const dictionary& dict,
57 const fvMesh& mesh
58)
59:
60 fv::option(sourceName, modelType, dict, mesh),
61 rhoName_(coeffs_.getOrDefault<word>("rho", "rho")),
62 Cmu_
63 (
64 dimensionedScalar::getOrAddToDict
65 (
66 "Cmu",
67 coeffs_,
68 0.09
69 )
70 ),
71 C_
72 (
73 dimensionedScalar::getOrAddToDict
74 (
75 "C",
76 coeffs_,
77 1.51
78 )
79 ),
80 lambda2_
81 (
82 dimensionedScalar::getOrAddToDict
83 (
84 "lambda2",
85 coeffs_,
86 0.1
87 )
88 ),
89 alpha_
90 (
91 dimensionedScalar::getOrAddToDict
92 (
93 "alpha",
94 coeffs_,
95 1.36
96 )
97 )
98{
100
101 // Note: incompressible only
102 const auto* turbPtr =
104 (
106 );
107
108 if (turbPtr)
109 {
110 const tmp<volScalarField>& tk = turbPtr->k();
111 fieldNames_[0] = tk().name();
112
113 const tmp<volScalarField>& tnut = turbPtr->nut();
114 fieldNames_[1] = tnut().name();
115
116 Log << " Applying model to " << fieldNames_[0]
117 << " and " << fieldNames_[1] << endl;
118 }
119 else
120 {
122 << "Unable to find incompressible turbulence model"
123 << exit(FatalError);
124 }
125
127}
128
129
130// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
131
133(
134 const volScalarField& rho,
135 fvMatrix<scalar>& eqn,
136 const label fieldi
137)
138{
139 // Not applicable to compressible cases
141}
142
143
145(
146 fvMatrix<scalar>& eqn,
147 const label fieldi
148)
149{
150 if (fieldi != 0)
151 {
152 return;
153 }
154
155 Log << this->name() << ": applying buoyancy production term to "
156 << eqn.psi().name() << endl;
157
158 // Buoyancy production in k eqn
159
160 const auto* turbPtr =
161 mesh_.findObject<incompressible::turbulenceModel>
162 (
164 );
165
166 if (!turbPtr)
167 {
169 << "Unable to find incompressible turbulence model"
170 << exit(FatalError);
171 }
172
173
174 tmp<volScalarField> tepsilon = turbPtr->epsilon();
175 const volScalarField& epsilon = tepsilon();
176 const volScalarField& k = eqn.psi();
177
178 // Note: using solver density field for incompressible multiphase cases
179 const auto& rho = mesh_.lookupObject<volScalarField>(rhoName_);
180
181 const auto& g = meshObjects::gravity::New(mesh_.time());
182
183 const dimensionedScalar eps0("eps0", epsilon.dimensions(), SMALL);
184
185 // Note: differing from reference by replacing nut/k by Cmu*k/epsilon
186 const volScalarField GbyK
187 (
188 "GbyK",
189 Cmu_*k/(epsilon + eps0)*alpha_*(g & fvc::grad(rho))/rho
190 );
191
192 return eqn -= fvm::SuSp(GbyK, k);
193}
194
195
197{
198 if (field.name() != fieldNames_[1])
199 {
200 return;
201 }
202
203 Log << this->name() << ": correcting " << field.name() << endl;
204
205 const auto* turbPtr =
206 mesh_.findObject<incompressible::turbulenceModel>
207 (
209 );
210
211 // nut correction
212 const auto& U = turbPtr->U();
213 tmp<volScalarField> tepsilon = turbPtr->epsilon();
214 const auto& epsilon = tepsilon();
215 tmp<volScalarField> tk = turbPtr->k();
216 const auto& k = tk();
217
219 const auto& gradU = tgradU();
220 const dimensionedScalar pSmall("pSmall", dimless/sqr(dimTime), SMALL);
221 const volScalarField pRat
222 (
223 magSqr(symm(gradU))/(magSqr(skew(gradU)) + pSmall)
224 );
225
226 const volScalarField epsilonTilde
227 (
228 max
229 (
230 epsilon,
231 lambda2_*C_*pRat*epsilon
232 )
233 );
234
235 field = Cmu_*sqr(k)/epsilonTilde;
236 field.correctBoundaryConditions();
237}
238
239
240// ************************************************************************* //
label k
#define Log
Definition: PDRblock.C:35
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
const uniformDimensionedVectorField & g
const word & name() const noexcept
Return the object name.
Definition: IOobjectI.H:65
Templated abstract base class for single-phase incompressible turbulence models.
void resize(const label len)
Adjust allocated size of list.
Definition: ListI.H:139
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
A special matrix type and solver, designed for finite volume solutions of scalar equations....
Definition: fvMatrix.H:121
const GeometricField< Type, fvPatchField, volMesh > & psi(const label i=0) const
Return psi.
Definition: fvMatrix.H:412
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
Applies corrections to the turbulent kinetic energy equation (i.e. k) and turbulent viscosity field (...
virtual void addSup(const volScalarField &rho, fvMatrix< scalar > &eqn, const label fieldi)
Add explicit contribution to compressible k equation.
Base abstract class for handling finite volume options (i.e. fvOption).
Definition: fvOption.H:127
const fvMesh & mesh_
Reference to the mesh database.
Definition: fvOption.H:139
wordList fieldNames_
Field names to apply source to - populated by derived models.
Definition: fvOption.H:148
void resetApplied()
Resize/reset applied flag list for all fieldNames_ entries.
Definition: fvOption.C:48
const Type * findObject(const word &name, const bool recursive=false) const
Return const pointer to the object of the given Type.
A class for managing temporary objects.
Definition: tmp.H:65
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from Foam::string.
Definition: word.H:68
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
U
Definition: pEqn.H:72
rDeltaTY field()
dynamicFvMesh & mesh
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:517
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
scalar epsilon
A special matrix type and solver, designed for finite volume solutions of scalar equations.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:54
zeroField SuSp(const Foam::zero, const GeometricField< Type, fvPatchField, volMesh > &)
A no-op source.
Namespace for OpenFOAM.
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
error FatalError
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:130
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
dimensionedTensor skew(const dimensionedTensor &dt)
labelList fv(nPoints)
dictionary dict