eddyViscosity.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) 2013-2023 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 
26 #include "eddyViscosity.H"
27 #include "fvcMeshPhi.H"
28 #include "fvmDiv.H"
29 #include "fvmSup.H"
30 
31 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32 
33 template<class BasicMomentumTransportModel>
35 (
36  const word& type,
37  const alphaField& alpha,
38  const rhoField& rho,
39  const volVectorField& U,
40  const surfaceScalarField& alphaRhoPhi,
41  const surfaceScalarField& phi,
42  const viscosity& viscosity
43 )
44 :
45  linearViscousStress<BasicMomentumTransportModel>
46  (
47  type,
48  alpha,
49  rho,
50  U,
51  alphaRhoPhi,
52  phi,
53  viscosity
54  ),
55 
56  nut_
57  (
58  IOobject
59  (
60  this->groupName("nut"),
61  this->runTime_.name(),
62  this->mesh_,
63  IOobject::MUST_READ,
64  IOobject::AUTO_WRITE
65  ),
66  this->mesh_
67  )
68 {}
69 
70 
71 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
72 
73 template<class BasicMomentumTransportModel>
75 {
77 }
78 
79 
80 template<class BasicMomentumTransportModel>
83 {
85  (
86  this->groupName("R"),
87  ((2.0/3.0)*I)*k() - (nut_)*dev(twoSymm(fvc::grad(this->U_)))
88  );
89 }
90 
91 
92 template<class BasicMomentumTransportModel>
94 {
95  correctNut();
96 }
97 
98 
99 template<class BasicMomentumTransportModel>
101 {
103 }
104 
105 
106 // ************************************************************************* //
label k
Generic GeometricField class.
static tmp< GeometricField< Type, PatchField, GeoMesh > > New(const word &name, const Internal &, const PtrList< PatchField< Type >> &, const HashPtrTable< Source > &=HashPtrTable< Source >())
Return a temporary field constructed from name,.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
BasicMomentumTransportModel::alphaField alphaField
Definition: eddyViscosity.H:70
virtual void validate()
Validate the turbulence fields after construction.
Definition: eddyViscosity.C:93
virtual bool read()=0
Re-read model coefficients if they have changed.
Definition: eddyViscosity.C:74
eddyViscosity(const word &modelName, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const viscosity &viscosity)
Construct from components.
Definition: eddyViscosity.C:35
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
virtual tmp< volSymmTensorField > sigma() const
Return the Reynolds stress tensor [m^2/s^2].
Definition: eddyViscosity.C:82
BasicMomentumTransportModel::rhoField rhoField
Definition: eddyViscosity.H:71
Linear viscous stress turbulence model base class.
A class for managing temporary objects.
Definition: tmp.H:55
Abstract base class for all fluid physical properties.
Definition: viscosity.H:50
A class for handling words, derived from string.
Definition: word.H:62
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
Calculate the matrix for the divergence of the given field and flux.
Calculate the matrix for implicit and explicit sources.
U
Definition: pEqn.H:72
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
void correct(const RdeltaTType &rDeltaT, const RhoType &rho, volScalarField &psi, const surfaceScalarField &phiCorr, const SpType &Sp, const SuType &Su)
void read(Istream &, label &, const dictionary &)
In-place read with dictionary lookup.
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
Definition: fvcGrad.C:46
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
static const Identity< scalar > I
Definition: Identity.H:93
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488