SmagorinskyZhang.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-2021 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 "SmagorinskyZhang.H"
27 #include "fvModels.H"
28 #include "fvConstraints.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace LESModels
35 {
36 
37 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
38 
39 template<class BasicMomentumTransportModel>
41 (
42  const alphaField& alpha,
43  const rhoField& rho,
44  const volVectorField& U,
45  const surfaceScalarField& alphaRhoPhi,
46  const surfaceScalarField& phi,
47  const viscosity& viscosity,
48  const word& type
49 )
50 :
52  (
53  alpha,
54  rho,
55  U,
56  alphaRhoPhi,
57  phi,
58  viscosity,
59  type
60  ),
61 
62  gasTurbulencePtr_(nullptr),
63 
64  Cmub_
65  (
67  (
68  "Cmub",
69  this->coeffDict_,
70  0.6
71  )
72  )
73 {
74  if (type == typeName)
75  {
76  this->printCoeffs(type);
77  }
78 }
79 
80 
81 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
82 
83 template<class BasicMomentumTransportModel>
85 {
87  {
88  Cmub_.readIfPresent(this->coeffDict());
89 
90  return true;
91  }
92  else
93  {
94  return false;
95  }
96 }
97 
98 
99 template<class BasicMomentumTransportModel>
102 {
103  if (!gasTurbulencePtr_)
104  {
105  const volVectorField& U = this->U_;
106 
107  const phaseModel& liquid =
108  refCast<const phaseModel>(this->properties());
109  const phaseSystem& fluid = liquid.fluid();
110  const phaseModel& gas = fluid.otherPhase(liquid);
111 
112  gasTurbulencePtr_ =
113  &U.db().lookupObject
114  <
116  >
117  (
119  (
120  momentumTransportModel::typeName,
121  gas.name()
122  )
123  );
124  }
125 
126  return *gasTurbulencePtr_;
127 }
128 
129 
130 template<class BasicMomentumTransportModel>
132 {
133  const phaseCompressibleMomentumTransportModel& gasTurbulence =
134  this->gasTurbulence();
135 
136  const phaseModel& liquid = refCast<const phaseModel>(this->properties());
137  const phaseSystem& fluid = liquid.fluid();
138  const phaseModel& gas = fluid.otherPhase(liquid);
139 
140  volScalarField k(this->k(fvc::grad(this->U_)));
141 
142  this->nut_ =
143  this->Ck_*sqrt(k)*this->delta()
144  + Cmub_*gas.d()*gasTurbulence.alpha()
145  *(mag(this->U_ - gasTurbulence.U()));
146 
147  this->nut_.correctBoundaryConditions();
148  fvConstraints::New(this->mesh_).constrain(this->nut_);
149 }
150 
151 
152 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
153 
154 } // End namespace LESModels
155 } // End namespace Foam
156 
157 // ************************************************************************* //
Generic thermophysical properties class for a liquid in which the functions and coefficients for each...
Definition: liquid.H:50
scalar delta
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
BasicMomentumTransportModel::rhoField rhoField
Templated abstract base class for multiphase compressible turbulence models.
The Smagorinsky SGS model including bubble-generated turbulence.
U
Definition: pEqn.H:72
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
const alphaField & alpha() const
Access function to phase fraction.
dimensionedScalar sqrt(const dimensionedScalar &ds)
const word & name() const
Definition: phaseModel.H:110
label k
Boltzmann constant.
Generic dimensioned Type class.
tmp< volScalarField > d() const
Return the Sauter-mean diameter.
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
The Smagorinsky SGS model.
Definition: Smagorinsky.H:86
Class to represent a system of phases and model interfacial transfers between them.
Definition: phaseSystem.H:68
A class for handling words, derived from string.
Definition: word.H:59
static word groupName(Name name, const word &group)
Info<< "Reading field p\"<< endl;volScalarField p(IOobject("p", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field U\"<< endl;volVectorField U(IOobject("U", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);pressureReference pressureReference(p, simple.dict());mesh.schemes().setFluxRequired(p.name());Info<< "Reading field pa\"<< endl;volScalarField pa(IOobject("pa", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field Ua\"<< endl;volVectorField Ua(IOobject("Ua", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);# 65 "/home/ubuntu/OpenFOAM-10/applications/solvers/incompressible/adjointShapeOptimisationFoam/createFields.H" 2label paRefCell=0;scalar paRefValue=0.0;setRefCell(pa, simple.dict(), paRefCell, paRefValue);mesh.schemes().setFluxRequired(pa.name());autoPtr< viscosityModel > viscosity(viscosityModel::New(mesh))
phaseSystem & fluid
Definition: createFields.H:11
phi
Definition: correctPhi.H:3
static autoPtr< dictionary > New(Istream &)
Construct top-level dictionary on freestore from Istream.
Definition: dictionaryIO.C:96
Abstract base class for all fluid physical properties.
Definition: viscosity.H:49
virtual bool read()
Read model coefficients if they have changed.
const volVectorField & U() const
Access function to velocity field.
BasicMomentumTransportModel::alphaField alphaField
SmagorinskyZhang(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const viscosity &viscosity, const word &type=typeName)
Construct from components.
const phaseSystem & fluid() const
Return the system to which this phase belongs.
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
dimensioned< scalar > mag(const dimensioned< Type > &)
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:312
virtual void correctNut()
Update the SGS eddy viscosity.
Single incompressible phase derived from the phase-fraction. Used as part of the multiPhaseMixture fo...
Definition: phaseModel.H:53
Namespace for OpenFOAM.