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 | Copyright (C) 2013-2015 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 "NicenoKEqn.H"
28 #include "twoPhaseSystem.H"
29 #include "dragModel.H"
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 namespace LESModels
36 {
37 
38 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
39 
40 template<class BasicTurbulenceModel>
41 NicenoKEqn<BasicTurbulenceModel>::NicenoKEqn
42 (
43  const alphaField& alpha,
44  const rhoField& rho,
45  const volVectorField& U,
46  const surfaceScalarField& alphaRhoPhi,
47  const surfaceScalarField& phi,
48  const transportModel& transport,
49  const word& propertiesName,
50  const word& type
51 )
52 :
54  (
55  alpha,
56  rho,
57  U,
58  alphaRhoPhi,
59  phi,
60  transport,
61  propertiesName,
62  type
63  ),
64 
65  gasTurbulencePtr_(NULL),
66 
67  alphaInversion_
68  (
70  (
71  "alphaInversion",
72  this->coeffDict_,
73  0.3
74  )
75  ),
76 
77  Cp_
78  (
80  (
81  "Cp",
82  this->coeffDict_,
83  this->Ck_.value()
84  )
85  ),
86 
87  Cmub_
88  (
90  (
91  "Cmub",
92  this->coeffDict_,
93  0.6
94  )
95  )
96 {
97  if (type == typeName)
98  {
99  // Cannot correct nut yet: construction of the phases is not complete
100  // correctNut();
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  else
121  {
122  return false;
123  }
124 }
125 
126 
127 template<class BasicTurbulenceModel>
129 <
130  typename BasicTurbulenceModel::transportModel
131 >&
133 {
134  if (!gasTurbulencePtr_)
135  {
136  const volVectorField& U = this->U_;
137 
138  const transportModel& liquid = this->transport();
139  const twoPhaseSystem& fluid =
140  refCast<const twoPhaseSystem>(liquid.fluid());
141  const transportModel& gas = fluid.otherPhase(liquid);
142 
143  gasTurbulencePtr_ =
144  &U.db()
146  (
148  (
150  gas.name()
151  )
152  );
153  }
154 
155  return *gasTurbulencePtr_;
156 }
157 
158 
159 template<class BasicTurbulenceModel>
161 {
163  this->gasTurbulence();
164 
165  this->nut_ =
166  this->Ck_*sqrt(this->k_)*this->delta()
167  + Cmub_*gasTurbulence.transport().d()*gasTurbulence.alpha()
168  *(mag(this->U_ - gasTurbulence.U()));
169 
170  this->nut_.correctBoundaryConditions();
171 }
172 
173 
174 template<class BasicTurbulenceModel>
176 {
178  this->gasTurbulence();
179 
180  const transportModel& liquid = this->transport();
181  const twoPhaseSystem& fluid =
182  refCast<const twoPhaseSystem>(liquid.fluid());
183  const transportModel& gas = fluid.otherPhase(liquid);
184 
185  volScalarField magUr(mag(this->U_ - gasTurbulence.U()));
186 
187  tmp<volScalarField> bubbleG
188  (
189  Cp_*sqr(magUr)*fluid.drag(gas).K()/liquid.rho()
190  );
191 
192  return bubbleG;
193 }
194 
195 
196 template<class BasicTurbulenceModel>
199 {
200  const volVectorField& U = this->U_;
201  const alphaField& alpha = this->alpha_;
202  const rhoField& rho = this->rho_;
203 
204  const turbulenceModel& gasTurbulence = this->gasTurbulence();
205 
206  return
207  (
208  max(alphaInversion_ - alpha, scalar(0))
209  *rho
210  *min
211  (
212  this->Ce_*sqrt(gasTurbulence.k())/this->delta(),
213  1.0/U.time().deltaT()
214  )
215  );
216 }
217 
218 
219 template<class BasicTurbulenceModel>
221 {
222  const alphaField& alpha = this->alpha_;
223  const rhoField& rho = this->rho_;
224 
226  this->gasTurbulence();
227 
228  const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
229 
230  return
231  alpha*rho*bubbleG()
232  + phaseTransferCoeff*gasTurbulence.k()
233  - fvm::Sp(phaseTransferCoeff, this->k_);
234 }
235 
236 
237 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
238 
239 } // End namespace LESModels
240 } // End namespace Foam
241 
242 // ************************************************************************* //
dimensionedScalar sqrt(const dimensionedScalar &ds)
virtual void correctNut()
Definition: NicenoKEqn.C:160
multiphaseSystem & fluid
Definition: createFields.H:10
One equation eddy-viscosity model.
Definition: kEqn.H:74
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual tmp< fvScalarMatrix > kSource() const
Definition: NicenoKEqn.C:220
BasicTurbulenceModel::rhoField rhoField
Definition: NicenoKEqn.H:123
A class for handling words, derived from string.
Definition: word.H:59
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
Namespace for OpenFOAM.
tmp< volScalarField > bubbleG() const
Definition: NicenoKEqn.C:175
static word groupName(Name name, const word &group)
Generic dimensioned Type class.
BasicTurbulenceModel::alphaField alphaField
Definition: NicenoKEqn.H:122
const phaseModel & otherPhase(const phaseModel &phase) const
Constant access the phase not given as an argument.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Class which solves the volume fraction equations for two phases.
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
One-equation SGS model for the continuous phase in a two-phase system including bubble-generated turb...
Definition: NicenoKEqn.H:75
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Templated abstract base class for multiphase compressible turbulence models.
scalar delta
Macros for easy insertion into run-time selection tables.
fileName::Type type(const fileName &)
Return the file type: DIRECTORY or FILE.
Definition: POSIX.C:589
surfaceScalarField & phi
Abstract base class for turbulence models (RAS, LES and laminar).
virtual tmp< volScalarField > K(const volScalarField &Ur) const =0
The dragfunction K used in the momentum eq.
tmp< volScalarField > phaseTransferCoeff() const
Definition: NicenoKEqn.C:198
BasicTurbulenceModel::transportModel transportModel
Definition: NicenoKEqn.H:124
const Time & time() const
Return time.
Definition: IOobject.C:251
static const word propertiesName
Default name of the turbulence properties dictionary.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
dimensionedScalar deltaT() const
Return time step.
Definition: TimeState.C:79
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:245
virtual bool read()
Read model coefficients if they have changed.
Definition: NicenoKEqn.C:110
const dragModel & drag(const phaseModel &phase) const
Return the drag model for the given phase.
A class for managing temporary objects.
Definition: PtrList.H:118
U
Definition: pEqn.H:82