kOmegaSST.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) 2007-2019 PCOpt/NTUA
9 Copyright (C) 2013-2019 FOSS GP
10 Copyright (C) 2019-2020 OpenCFD Ltd.
11-------------------------------------------------------------------------------
12License
13 This file is part of OpenFOAM.
14
15 OpenFOAM is free software: you can redistribute it and/or modify it
16 under the terms of the GNU General Public License as published by
17 the Free Software Foundation, either version 3 of the License, or
18 (at your option) any later version.
19
20 OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
21 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
22 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
23 for more details.
24
25 You should have received a copy of the GNU General Public License
26 along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
27
28\*---------------------------------------------------------------------------*/
29
30#include "kOmegaSST.H"
31#include "wallDist.H"
33
34// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35
36namespace Foam
37{
38namespace incompressible
39{
40namespace RASVariables
41{
42
43// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
44
47
48
49// * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
50
52{
55 {
56 GMean_.reset
57 (
59 (
61 (
62 "GMean",
64 mesh_,
67 ),
68 mesh_,
70 )
71 );
72 }
73}
74
75
76// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
77
79(
80 const fvMesh& mesh,
81 const solverControl& SolverControl
82)
83:
84 RASModelVariables(mesh, SolverControl)
85{
86 TMVar1BaseName_ = "k";
87 TMVar2BaseName_ = "omega";
88
92
93 distPtr_.ref(const_cast<volScalarField&>(wallDist::New(mesh_).y()));
94
97}
98
99
101{
103 (
105 (
107 TMVar2().internalField().group()
108 )
109 );
110 // Recompute G and modify values next to the walls
111 // Ideally, grad(U) should be cached to avoid the overhead
112 const volVectorField& U = turbModel.U();
115 (
116 this->type() + ":GbyNu",
117 (tgradU() && dev(twoSymm(tgradU())))
118 );
119 auto tG =
121 (
122 turbModel.GName(),
123 nutRefInst()*GbyNu0
124 );
125 // Use correctBoundaryConditions instead of updateCoeffs to avoid
126 // messing with updateCoeffs in the next iteration of omegaEqn
128
129 return tG;
130}
131
132
133// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
134
136{
138 {
140 << "Using GMean" << endl;
142 }
144 << "Using instantaneous G" << endl;
145 return computeG();
146}
147
148
150{
153 {
154 const label iAverageIter = solverControl_.averageIter();
155 scalar avIter(iAverageIter);
156 scalar oneOverItP1 = 1./(avIter + 1);
157 scalar mult = avIter*oneOverItP1;
158 GMean_() = GMean_()*mult + computeG()*oneOverItP1;
159 }
160}
161
162
164(
165 const incompressible::turbulenceModel& turbulence
166)
167{
168 // The presence of G is required to update the boundary value of omega
169 const volVectorField& U = turbulence.U();
170 const volScalarField S2(2*magSqr(symm(fvc::grad(U))));
171
172 volScalarField G(turbulence.GName(), nutRef() * S2);
174}
175
176
177// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
178
179} // End namespace RASVariables
180} // End namespace incompressible
181} // End namespace Foam
182
183// ************************************************************************* //
scalar y
Macros for easy insertion into run-time selection tables.
#define addToRunTimeSelectionTable(baseType, thisType, argNames)
Add to construction table with typeName as the key.
void correctBoundaryConditions()
Correct boundary field.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
Templated abstract base class for single-phase incompressible turbulence models.
static autoPtr< Time > New()
Construct (dummy) Time - no functionObjects or libraries.
Definition: Time.C:717
static word timeName(const scalar t, const int precision=precision_)
Definition: Time.C:780
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:126
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:91
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:290
Abstract base class for objective functions. No point in making this runTime selectable since its chi...
const volScalarField & TMVar2() const
const volScalarField & nutRefInst() const
virtual void computeMeanFields()
Compute mean fields on the fly.
const volScalarField & nutRef() const
const volScalarField & TMVar2Inst() const
tmp< volScalarField::Internal > computeG()
Definition: kOmegaSST.C:100
virtual tmp< volScalarField::Internal > G()
Return the turbulence production term.
Definition: kOmegaSST.C:135
virtual void computeMeanFields()
Compute mean fields on the fly.
Definition: kOmegaSST.C:149
autoPtr< volScalarField::Internal > GMean_
Average of the production term.
Definition: kOmegaSST.H:67
Type & lookupObjectRef(const word &name, const bool recursive=false) const
const Type & lookupObject(const word &name, const bool recursive=false) const
static const word propertiesName
Default name of the phase properties dictionary.
Definition: phaseSystem.H:290
Base class for solver control classes.
Definition: solverControl.H:52
bool doAverageIter() const
label & averageIter()
Return average iteration index reference.
bool useAveragedFields() const
bool average() const
Whether averaging is enabled or not.
A class for managing temporary objects.
Definition: tmp.H:65
const volVectorField & U() const
Access function to velocity field.
word GName() const
Helper function to return the name of the turbulence G field.
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
U
Definition: pEqn.H:72
dynamicFvMesh & mesh
compressible::turbulenceModel & turbulence
#define DebugInfo
Report an information message using Foam::Info.
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh > > grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:54
Namespace for OpenFOAM.
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
dimensionedScalar pow3(const dimensionedScalar &ds)
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:53
const dimensionSet dimArea(sqr(dimLength))
Definition: dimensionSets.H:59
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:598
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:372
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)