NicenoKEqn.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) 2013-2016 OpenFOAM Foundation
9  Copyright (C) 2019-2020 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
29 #include "NicenoKEqn.H"
30 #include "fvOptions.H"
31 #include "twoPhaseSystem.H"
32 #include "dragModel.H"
33 
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 namespace LESModels
39 {
40 
41 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
42 
43 template<class BasicTurbulenceModel>
44 NicenoKEqn<BasicTurbulenceModel>::NicenoKEqn
45 (
46  const alphaField& alpha,
47  const rhoField& rho,
48  const volVectorField& U,
49  const surfaceScalarField& alphaRhoPhi,
50  const surfaceScalarField& phi,
51  const transportModel& transport,
52  const word& propertiesName,
53  const word& type
54 )
55 :
57  (
58  alpha,
59  rho,
60  U,
61  alphaRhoPhi,
62  phi,
63  transport,
64  propertiesName,
65  type
66  ),
67 
68  gasTurbulencePtr_(nullptr),
69 
70  alphaInversion_
71  (
73  (
74  "alphaInversion",
75  this->coeffDict_,
76  0.3
77  )
78  ),
79 
80  Cp_
81  (
83  (
84  "Cp",
85  this->coeffDict_,
86  this->Ck_.value()
87  )
88  ),
89 
90  Cmub_
91  (
93  (
94  "Cmub",
95  this->coeffDict_,
96  0.6
97  )
98  )
99 {
100  if (type == typeName)
101  {
102  this->printCoeffs(type);
103  }
104 }
105 
106 
107 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
108 
109 template<class BasicTurbulenceModel>
111 {
113  {
114  alphaInversion_.readIfPresent(this->coeffDict());
115  Cp_.readIfPresent(this->coeffDict());
116  Cmub_.readIfPresent(this->coeffDict());
117 
118  return true;
119  }
120 
121  return false;
122 }
123 
124 
125 template<class BasicTurbulenceModel>
127 <
128  typename BasicTurbulenceModel::transportModel
129 >&
131 {
132  if (!gasTurbulencePtr_)
133  {
134  const volVectorField& U = this->U_;
135 
136  const transportModel& liquid = this->transport();
137  const twoPhaseSystem& fluid =
138  refCast<const twoPhaseSystem>(liquid.fluid());
139  const transportModel& gas = fluid.otherPhase(liquid);
140 
141  gasTurbulencePtr_ =
142  &U.db()
144  (
146  (
148  gas.name()
149  )
150  );
151  }
152 
153  return *gasTurbulencePtr_;
154 }
155 
156 
157 template<class BasicTurbulenceModel>
159 {
161  this->gasTurbulence();
162 
163  this->nut_ =
164  this->Ck_*sqrt(this->k_)*this->delta()
165  + Cmub_*gasTurbulence.transport().d()*gasTurbulence.alpha()
166  *(mag(this->U_ - gasTurbulence.U()));
167 
168  this->nut_.correctBoundaryConditions();
169  fv::options::New(this->mesh_).correct(this->nut_);
170 
171  BasicTurbulenceModel::correctNut();
172 }
173 
174 
175 template<class BasicTurbulenceModel>
177 {
179  this->gasTurbulence();
180 
181  const transportModel& liquid = this->transport();
182  const twoPhaseSystem& fluid =
183  refCast<const twoPhaseSystem>(liquid.fluid());
184  const transportModel& gas = fluid.otherPhase(liquid);
185 
186  const dragModel& drag = fluid.lookupSubModel<dragModel>(gas, liquid);
187 
188  volScalarField magUr(mag(this->U_ - gasTurbulence.U()));
189 
190  tmp<volScalarField> bubbleG
191  (
192  Cp_*sqr(magUr)*drag.K()/liquid.rho()
193  );
194 
195  return bubbleG;
196 }
197 
198 
199 template<class BasicTurbulenceModel>
202 {
203  const volVectorField& U = this->U_;
204  const alphaField& alpha = this->alpha_;
205  const rhoField& rho = this->rho_;
206 
207  const turbulenceModel& gasTurbulence = this->gasTurbulence();
208 
209  return
210  (
211  max(alphaInversion_ - alpha, scalar(0))
212  *rho
213  *min
214  (
215  this->Ce_*sqrt(gasTurbulence.k())/this->delta(),
216  1.0/U.time().deltaT()
217  )
218  );
219 }
220 
221 
222 template<class BasicTurbulenceModel>
224 {
225  const alphaField& alpha = this->alpha_;
226  const rhoField& rho = this->rho_;
227 
229  this->gasTurbulence();
230 
231  const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
232 
233  return
234  alpha*rho*bubbleG()
235  + phaseTransferCoeff*gasTurbulence.k()
236  - fvm::Sp(phaseTransferCoeff, this->k_);
237 }
238 
239 
240 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
241 
242 } // End namespace LESModels
243 } // End namespace Foam
244 
245 // ************************************************************************* //
drag
Info<< "Reading strained laminar flame speed field Su\n"<< endl;volScalarField Su(IOobject("Su", runTime.timeName(), mesh, IOobject::MUST_READ, IOobject::AUTO_WRITE), mesh);Info<< "Reading field betav\n"<< endl;volScalarField betav(IOobject("betav", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Lobs\n"<< endl;volScalarField Lobs(IOobject("Lobs", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field CT\n"<< endl;volSymmTensorField CT(IOobject("CT", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field Nv\n"<< endl;volScalarField Nv(IOobject("Nv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);Info<< "Reading field nsv\n"<< endl;volSymmTensorField nsv(IOobject("nsv", mesh.facesInstance(), mesh, IOobject::MUST_READ, IOobject::NO_WRITE), mesh);IOdictionary PDRProperties(IOobject("PDRProperties", runTime.constant(), mesh, IOobject::MUST_READ_IF_MODIFIED, IOobject::NO_WRITE));autoPtr< PDRDragModel > drag
Definition: createFields.H:165
fvOptions.H
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::twoPhaseSystem
Class which solves the volume fraction equations for two phases.
Definition: twoPhaseSystem.H:53
Foam::fv::optionList::correct
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
Definition: fvOptionListTemplates.C:355
Foam::PhaseCompressibleTurbulenceModel
Templated abstract base class for multiphase compressible turbulence models.
Definition: phaseCompressibleTurbulenceModelFwd.H:44
Foam::tmp
A class for managing temporary objects.
Definition: PtrList.H:61
Foam::LESModels::NicenoKEqn::phaseTransferCoeff
tmp< volScalarField > phaseTransferCoeff() const
Definition: NicenoKEqn.C:201
Foam::fv::options::New
static options & New(const fvMesh &mesh)
Construct fvOptions and register to database if not present.
Definition: fvOptions.C:103
Foam::constant::atomic::alpha
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Definition: readThermalProperties.H:212
Foam::turbulenceModel::propertiesName
static const word propertiesName
Default name of the turbulence properties dictionary.
Definition: turbulenceModel.H:100
Foam::LESModels::NicenoKEqn::rhoField
BasicTurbulenceModel::rhoField rhoField
Definition: NicenoKEqn.H:128
Foam::liquid
Generic thermophysical properties class for a liquid in which the functions and coefficients for each...
Definition: liquid.H:54
fluid
twoPhaseSystem & fluid
Definition: setRegionFluidFields.H:3
rho
rho
Definition: readInitialConditions.H:88
Foam::min
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
Foam::LESModels::kEqn
One equation eddy-viscosity model.
Definition: kEqn.H:77
NicenoKEqn.H
Foam::LESModels::NicenoKEqn::kSource
virtual tmp< fvScalarMatrix > kSource() const
Definition: NicenoKEqn.C:223
Foam::LESModels::NicenoKEqn::alphaField
BasicTurbulenceModel::alphaField alphaField
Definition: NicenoKEqn.H:127
Foam::fvm::Sp
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
Foam::LESModels::NicenoKEqn
One-equation SGS model for the continuous phase in a two-phase system including bubble-generated turb...
Definition: NicenoKEqn.H:78
delta
scalar delta
Definition: LISASMDCalcMethod2.H:8
phi
surfaceScalarField & phi
Definition: setRegionFluidFields.H:8
Foam::LESModels::NicenoKEqn::read
virtual bool read()
Read model coefficients if they have changed.
Definition: NicenoKEqn.C:110
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
Foam::LESModels::NicenoKEqn::bubbleG
tmp< volScalarField > bubbleG() const
Definition: NicenoKEqn.C:176
Foam::turbulenceModel
Abstract base class for turbulence models (RAS, LES and laminar).
Definition: turbulenceModel.H:63
Foam::dimensioned
Generic dimensioned Type class.
Definition: dimensionedScalarFwd.H:42
Foam::turbulenceModel::k
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::LESModels::NicenoKEqn::correctNut
virtual void correctNut()
Definition: NicenoKEqn.C:158
Foam::transportModel
Base-class for all transport models used by the incompressible turbulence models.
Definition: transportModel.H:53
Foam::liquid::rho
scalar rho(scalar p, scalar T) const
Liquid density [kg/m^3].
Definition: liquidI.H:28
U
U
Definition: pEqn.H:72
Foam::sqr
dimensionedSymmTensor sqr(const dimensionedVector &dv)
Definition: dimensionedSymmTensor.C:51
Foam::LESModels::NicenoKEqn::transportModel
BasicTurbulenceModel::transportModel transportModel
Definition: NicenoKEqn.H:129
Foam::sqrt
dimensionedScalar sqrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:144
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::type
fileName::Type type(const fileName &name, const bool followLink=true)
Return the file type: DIRECTORY or FILE, normally following symbolic links.
Definition: MSwindows.C:590
Foam::dragModel
Definition: dragModel.H:51
Foam::IOobject::groupName
static word groupName(StringType base, const word &group)
Create dot-delimited name.group string.
Foam::GeometricField< vector, fvPatchField, volMesh >