kEpsilonPhitF.H 8.58 KB
Newer Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
/*---------------------------------------------------------------------------*\
  =========                 |
  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   \\    /   O peration     |
    \\  /    A nd           | www.openfoam.com
     \\/     M anipulation  |
-------------------------------------------------------------------------------
    Copyright (C) 2019 OpenCFD Ltd.
-------------------------------------------------------------------------------
License
    This file is part of OpenFOAM.

    OpenFOAM is free software: you can redistribute it and/or modify it
    under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
    for more details.

    You should have received a copy of the GNU General Public License
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.

Class
    Foam::RASModels::kEpsilonPhitF

Group
    grpRASTurbulence

Description
    The k-epsilon-phit-f turbulence closure model for incompressible and
    compressible flows.

    The model is a three-transport-equation linear-eddy-viscosity turbulence
    closure model alongside an elliptic relaxation equation:
      - Turbulent kinetic energy, \c k,
      - Turbulent kinetic energy dissipation rate, \c epsilon,
      - Normalised wall-normal fluctuating velocity scale, \c phit,
      - Elliptic relaxation factor, \c f.

    Reference:
    \verbatim
        Standard model (Tag:LUU):
            Laurence, D. R., Uribe, J. C., & Utyuzhnikov, S. V. (2005).
47
            A robust formulation of the v2-f model.
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
            Flow, Turbulence and Combustion, 73(3-4), 169–185.
            DOI:10.1007/s10494-005-1974-8
    \endverbatim

    The default model coefficients are (LUU:Eqs. 19-20):
    \verbatim
        kEpsilonPhitFCoeffs
        {
            Cmu         0.22,    // Turbulent viscosity constant
            Ceps1a      1.4,     // Model constant for epsilon
            Ceps1b      1.0,     // Model constant for epsilon
            Ceps1c      0.05,    // Model constant for epsilon
            Ceps2       1.9,     // Model constant for epsilon
            Cf1         1.4,     // Model constant for f
            Cf2         0.3,     // Model constant for f
            CL          0.25,    // Model constant for L
            Ceta        110.0,   // Model constant for L
            CT          6.0,     // Model constant for T
            sigmaK      1.0,     // Turbulent Prandtl number for k
            sigmaEps    1.3,     // Turbulent Prandtl number for epsilon
            sigmaPhit   1.0,     // Turbulent Prandtl number for phit = sigmaK
        }
    \endverbatim

Note
    The name of the original variable replacing 'v2' is 'phi' (LUU:Eq. 14).
    However, the name 'phi' preexisted in OpenFOAM; therefore, this name was
    replaced by 'phit' herein.

SourceFiles
    kEpsilonPhitF.C

SeeAlso
    v2f.H

\*---------------------------------------------------------------------------*/

#ifndef kEpsilonPhitF_H
#define kEpsilonPhitF_H

#include "kEpsilonPhitFBase.H"
#include "RASModel.H"
#include "eddyViscosity.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

namespace Foam
{
namespace RASModels
{

/*---------------------------------------------------------------------------*\
                        Class kEpsilonPhitF Declaration
\*---------------------------------------------------------------------------*/

template<class BasicTurbulenceModel>
class kEpsilonPhitF
:
    public eddyViscosity<RASModel<BasicTurbulenceModel>>,
    public kEpsilonPhitFBase
{
    // Private Member Functions

        //- No copy construct
        kEpsilonPhitF(const kEpsilonPhitF&) = delete;

        //- No copy assignment
        void operator=(const kEpsilonPhitF&) = delete;


protected:

    // Protected Data

        // Model coefficients

            dimensionedScalar Cmu_;
            dimensionedScalar Ceps1a_;
            dimensionedScalar Ceps1b_;
            dimensionedScalar Ceps1c_;
            dimensionedScalar Ceps2_;
            dimensionedScalar Cf1_;
            dimensionedScalar Cf2_;
            dimensionedScalar CL_;
            dimensionedScalar Ceta_;
            dimensionedScalar CT_;
            dimensionedScalar sigmaK_;
            dimensionedScalar sigmaEps_;
            dimensionedScalar sigmaPhit_;


        // Fields

            //- Turbulent kinetic energy [m2/s2]
            volScalarField k_;

            //- Turbulent kinetic energy dissipation rate [m2/s3]
            volScalarField epsilon_;

            //- Normalised wall-normal fluctuating velocity scale [-]
            volScalarField phit_;

            //- Elliptic relaxation factor [1/s]
            volScalarField f_;

            //- Turbulent time scale [s]
            volScalarField T_;


        // Bounding values

            dimensionedScalar phitMin_;
            dimensionedScalar fMin_;
            dimensionedScalar TMin_;
            dimensionedScalar L2Min_;


    // Protected Member Functions

        //- Update nut with the latest available k, phit, and T
        virtual void correctNut();

        //- Return the turbulent time scale, T
171
        tmp<volScalarField> Ts() const;
172
173

        //- Return the turbulent length scale, L
174
        tmp<volScalarField> Ls() const;
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295


public:

    typedef typename BasicTurbulenceModel::alphaField alphaField;
    typedef typename BasicTurbulenceModel::rhoField rhoField;
    typedef typename BasicTurbulenceModel::transportModel transportModel;


    //- Runtime type information
    TypeName("kEpsilonPhitF");


    // Constructors

        //- Construct from components
        kEpsilonPhitF
        (
            const alphaField& alpha,
            const rhoField& rho,
            const volVectorField& U,
            const surfaceScalarField& alphaRhoPhi,
            const surfaceScalarField& phi,
            const transportModel& transport,
            const word& propertiesName = turbulenceModel::propertiesName,
            const word& type = typeName
        );


    //- Destructor
    virtual ~kEpsilonPhitF() = default;


    // Member Functions

        //- Re-read model coefficients if they have changed
        virtual bool read();

        //- Return the effective diffusivity for k
        tmp<volScalarField> DkEff() const
        {
            return tmp<volScalarField>
            (
                new volScalarField
                (
                    "DkEff",
                    this->nut_/sigmaK_ + this->nu()
                )
            );
        }

        //- Return the effective diffusivity for epsilon
        tmp<volScalarField> DepsilonEff() const
        {
            return tmp<volScalarField>
            (
                new volScalarField
                (
                    "DepsilonEff",
                    this->nut_/sigmaEps_ + this->nu()
                )
            );
        }

        //- Return the effective diffusivity for phit
        tmp<volScalarField> DphitEff() const
        {
            return tmp<volScalarField>
            (
                new volScalarField
                (
                    "DphitEff",
                    this->nut_/sigmaPhit_ + this->nu()
                )
            );
        }

        //- Return the turbulent kinetic energy field
        virtual tmp<volScalarField> k() const
        {
            return k_;
        }

        //- Return the turbulent kinetic energy dissipation rate field
        virtual tmp<volScalarField> epsilon() const
        {
            return epsilon_;
        }

        //- Return the normalised wall-normal fluctuating velocity scale field
        virtual tmp<volScalarField> phit() const
        {
            return phit_;
        }

        //- Return the elliptic relaxation factor field
        virtual tmp<volScalarField> f() const
        {
            return f_;
        }

        //- Solve the transport equations and correct the turbulent viscosity
        virtual void correct();
};


// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

} // End namespace RASModels
} // End namespace Foam

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#ifdef NoRepository
    #include "kEpsilonPhitF.C"
#endif

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

#endif

296
// ************************************************************************* //