SpalartAllmaras.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 "SpalartAllmaras.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// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
49
51(
52 const fvMesh& mesh,
53 const solverControl& SolverControl
54)
55:
56 RASModelVariables(mesh, SolverControl)
57{
58 TMVar1BaseName_ = "nuTilda";
59
61
62 TMVar2Ptr_.reset
63 (
65 (
67 (
68 "dummySpalartAllmarasVar2",
69 mesh.time().timeName(),
70 mesh,
73 ),
74 mesh,
76 )
77 );
78
80
81 // The wall dist name can vary depending on how wallDist was
82 // constructed. Grab the field directly from wallDist
83
84 distPtr_.ref
85 (
86 const_cast<volScalarField&>(wallDist::New(mesh_).y())
87 );
88
91}
92
93
94// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
95
97(
99) const
100{
101 auto tnutJacobian = tmp<volScalarField>::New
102 (
104 (
105 "nutJacobianVar1",
106 mesh_.time().timeName(),
107 mesh_,
110 ),
111 mesh_,
113 );
114 auto& nutJacobian = tnutJacobian.ref();
115
117 const volScalarField& nuTilda = TMVar1();
118
119 volScalarField chi(nuTilda/nu);
120 volScalarField chi3(pow3(chi));
121
122 const scalar Cv13 = pow3(7.1);
123 volScalarField fv1(chi3/(chi3 + Cv13));
124 volScalarField fv1ByChi2Sqr(sqr(chi/(chi3 + Cv13)));
125 volScalarField Cdfv1(3.0*Cv13*fv1ByChi2Sqr);
126
127 nutJacobian = Cdfv1*chi + fv1;
128
129 return tnutJacobian;
130}
131
132
133// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
134
135} // End namespace RASVariables
136} // End namespace incompressible
137} // End namespace Foam
138
139// ************************************************************************* //
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.
Defines the attributes of an object for which implicit objectRegistry management is supported,...
Definition: IOobject.H:170
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 & TMVar1() const
Return references to turbulence fields.
virtual tmp< volScalarField > nutJacobianVar1(const singlePhaseTransportModel &laminarTransport) const
return nut Jacobian wrt the TM vars
Type & lookupObjectRef(const word &name, const bool recursive=false) const
A simple single-phase transport model based on viscosityModel.
virtual tmp< volScalarField > nu() const
Return the laminar viscosity.
Base class for solver control classes.
Definition: solverControl.H:52
A class for managing temporary objects.
Definition: tmp.H:65
#define defineTypeNameAndDebug(Type, DebugSwitch)
Define the typeName and debug information.
Definition: className.H:121
dynamicFvMesh & mesh
Namespace for OpenFOAM.
const dimensionSet dimless
Dimensionless.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow3(const dimensionedScalar &ds)
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
volScalarField & nu
singlePhaseTransportModel laminarTransport(U, phi)