continuousGasKEqn.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-2018 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 "continuousGasKEqn.H"
27 #include "twoPhaseSystem.H"
28 
29 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
30 
31 namespace Foam
32 {
33 namespace LESModels
34 {
35 
36 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
37 
38 template<class BasicTurbulenceModel>
40 (
41  const alphaField& alpha,
42  const rhoField& rho,
43  const volVectorField& U,
44  const surfaceScalarField& alphaRhoPhi,
45  const surfaceScalarField& phi,
46  const transportModel& transport,
47  const word& propertiesName,
48  const word& type
49 )
50 :
52  (
53  alpha,
54  rho,
55  U,
56  alphaRhoPhi,
57  phi,
58  transport,
59  propertiesName,
60  type
61  ),
62 
63  liquidTurbulencePtr_(nullptr),
64 
65  alphaInversion_
66  (
68  (
69  "alphaInversion",
70  this->coeffDict_,
71  0.7
72  )
73  )
74 {
75  if (type == typeName)
76  {
77  this->printCoeffs(type);
78  }
79 }
80 
81 
82 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
83 
84 template<class BasicTurbulenceModel>
86 {
88  {
89  alphaInversion_.readIfPresent(this->coeffDict());
90 
91  return true;
92  }
93  else
94  {
95  return false;
96  }
97 }
98 
99 
100 template<class BasicTurbulenceModel>
101 const turbulenceModel&
103 {
104  if (!liquidTurbulencePtr_)
105  {
106  const volVectorField& U = this->U_;
107 
108  const transportModel& gas = this->transport();
109  const twoPhaseSystem& fluid =
110  refCast<const twoPhaseSystem>(gas.fluid());
111  const transportModel& liquid = fluid.otherPhase(gas);
112 
113  liquidTurbulencePtr_ =
115  (
117  (
119  liquid.name()
120  )
121  );
122  }
123 
124  return *liquidTurbulencePtr_;
125 }
126 
127 
128 template<class BasicTurbulenceModel>
131 {
132  const volVectorField& U = this->U_;
133  const alphaField& alpha = this->alpha_;
134  const rhoField& rho = this->rho_;
135 
136  const turbulenceModel& liquidTurbulence = this->liquidTurbulence();
137 
138  return
139  (
140  max(alphaInversion_ - alpha, scalar(0))
141  *rho
142  *min
143  (
144  this->Ce_*sqrt(liquidTurbulence.k())/this->delta(),
145  1.0/U.time().deltaT()
146  )
147  );
148 }
149 
150 
151 template<class BasicTurbulenceModel>
154 {
155  const turbulenceModel& liquidTurbulence = this->liquidTurbulence();
156  const volScalarField phaseTransferCoeff(this->phaseTransferCoeff());
157 
158  return
159  phaseTransferCoeff*liquidTurbulence.k()
160  - fvm::Sp(phaseTransferCoeff, this->k_);
161 }
162 
163 
164 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
165 
166 } // End namespace LESModels
167 } // End namespace Foam
168 
169 // ************************************************************************* //
scalar delta
surfaceScalarField & phi
multiphaseSystem & fluid
Definition: createFields.H:11
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedScalar sqrt(const dimensionedScalar &ds)
Generic dimensioned Type class.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
Abstract base class for turbulence models (RAS, LES and laminar).
One equation eddy-viscosity model.
Definition: kEqn.H:71
const turbulenceModel & liquidTurbulence() const
Return the turbulence model for the liquid phase.
Class which solves the volume fraction equations for two phases.
virtual tmp< volScalarField > k() const =0
Return the turbulence kinetic energy.
static const word propertiesName
Default name of the turbulence properties dictionary.
A class for handling words, derived from string.
Definition: word.H:59
static word groupName(Name name, const word &group)
BasicTurbulenceModel::alphaField alphaField
const phaseModel & otherPhase(const phaseModel &phase) const
Return the phase not given as an argument.
virtual tmp< fvScalarMatrix > kSource() const
tmp< volScalarField > phaseTransferCoeff() const
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
virtual bool read()
Read model coefficients if they have changed.
U
Definition: pEqn.H:72
BasicTurbulenceModel::transportModel transportModel
const Time & time() const
Return time.
Definition: IOobject.C:360
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
continuousGasKEqn(const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName=turbulenceModel::propertiesName, const word &type=typeName)
Construct from components.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
A class for managing temporary objects.
Definition: PtrList.H:53
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:354
BasicTurbulenceModel::rhoField rhoField
dimensionedScalar deltaT() const
Return time step.
Definition: TimeStateI.H:53
Namespace for OpenFOAM.