kOmega.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) 2011-2019 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 "kOmega.H"
27 #include "fvOptions.H"
28 #include "bound.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace RASModels
35 {
36 
37 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
38 
39 template<class BasicTurbulenceModel>
41 {
42  this->nut_ = k_/omega_;
43  this->nut_.correctBoundaryConditions();
44  fv::options::New(this->mesh_).correct(this->nut_);
45 
46  BasicTurbulenceModel::correctNut();
47 }
48 
49 
50 template<class BasicTurbulenceModel>
52 {
53  return tmp<fvScalarMatrix>
54  (
55  new fvScalarMatrix
56  (
57  k_,
58  dimVolume*this->rho_.dimensions()*k_.dimensions()
59  /dimTime
60  )
61  );
62 }
63 
64 
65 template<class BasicTurbulenceModel>
67 {
68  return tmp<fvScalarMatrix>
69  (
70  new fvScalarMatrix
71  (
72  omega_,
73  dimVolume*this->rho_.dimensions()*omega_.dimensions()/dimTime
74  )
75  );
76 }
77 
78 
79 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
80 
81 template<class BasicTurbulenceModel>
83 (
84  const alphaField& alpha,
85  const rhoField& rho,
86  const volVectorField& U,
87  const surfaceScalarField& alphaRhoPhi,
88  const surfaceScalarField& phi,
89  const transportModel& transport,
90  const word& propertiesName,
91  const word& type
92 )
93 :
95  (
96  type,
97  alpha,
98  rho,
99  U,
100  alphaRhoPhi,
101  phi,
102  transport,
103  propertiesName
104  ),
105 
106  Cmu_
107  (
109  (
110  "betaStar",
111  this->coeffDict_,
112  0.09
113  )
114  ),
115  beta_
116  (
118  (
119  "beta",
120  this->coeffDict_,
121  0.072
122  )
123  ),
124  gamma_
125  (
127  (
128  "gamma",
129  this->coeffDict_,
130  0.52
131  )
132  ),
133  alphaK_
134  (
136  (
137  "alphaK",
138  this->coeffDict_,
139  0.5
140  )
141  ),
142  alphaOmega_
143  (
145  (
146  "alphaOmega",
147  this->coeffDict_,
148  0.5
149  )
150  ),
151 
152  k_
153  (
154  IOobject
155  (
156  IOobject::groupName("k", alphaRhoPhi.group()),
157  this->runTime_.timeName(),
158  this->mesh_,
161  ),
162  this->mesh_
163  ),
164  omega_
165  (
166  IOobject
167  (
168  IOobject::groupName("omega", alphaRhoPhi.group()),
169  this->runTime_.timeName(),
170  this->mesh_,
173  ),
174  this->mesh_
175  )
176 {
177  bound(k_, this->kMin_);
178  bound(omega_, this->omegaMin_);
179 
180  if (type == typeName)
181  {
182  this->printCoeffs(type);
183  }
184 }
185 
186 
187 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
188 
189 template<class BasicTurbulenceModel>
191 {
193  {
194  Cmu_.readIfPresent(this->coeffDict());
195  beta_.readIfPresent(this->coeffDict());
196  gamma_.readIfPresent(this->coeffDict());
197  alphaK_.readIfPresent(this->coeffDict());
198  alphaOmega_.readIfPresent(this->coeffDict());
199 
200  return true;
201  }
202  else
203  {
204  return false;
205  }
206 }
207 
208 
209 template<class BasicTurbulenceModel>
211 {
212  if (!this->turbulence_)
213  {
214  return;
215  }
216 
217  // Local references
218  const alphaField& alpha = this->alpha_;
219  const rhoField& rho = this->rho_;
220  const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi_;
221  const volVectorField& U = this->U_;
222  volScalarField& nut = this->nut_;
223  fv::options& fvOptions(fv::options::New(this->mesh_));
224 
226 
228  (
229  fvc::div(fvc::absolute(this->phi(), U))().v()
230  );
231 
232  tmp<volTensorField> tgradU = fvc::grad(U);
234  (
235  this->GName(),
236  nut.v()*(dev(twoSymm(tgradU().v())) && tgradU().v())
237  );
238  tgradU.clear();
239 
240  // Update omega and G at the wall
241  omega_.boundaryFieldRef().updateCoeffs();
242 
243  // Turbulence specific dissipation rate equation
244  tmp<fvScalarMatrix> omegaEqn
245  (
246  fvm::ddt(alpha, rho, omega_)
247  + fvm::div(alphaRhoPhi, omega_)
248  - fvm::laplacian(alpha*rho*DomegaEff(), omega_)
249  ==
250  gamma_*alpha()*rho()*G*omega_()/k_()
251  - fvm::SuSp(((2.0/3.0)*gamma_)*alpha()*rho()*divU, omega_)
252  - fvm::Sp(beta_*alpha()*rho()*omega_(), omega_)
253  + omegaSource()
254  + fvOptions(alpha, rho, omega_)
255  );
256 
257  omegaEqn.ref().relax();
258  fvOptions.constrain(omegaEqn.ref());
259  omegaEqn.ref().boundaryManipulate(omega_.boundaryFieldRef());
260  solve(omegaEqn);
261  fvOptions.correct(omega_);
262  bound(omega_, this->omegaMin_);
263 
264 
265  // Turbulent kinetic energy equation
267  (
268  fvm::ddt(alpha, rho, k_)
269  + fvm::div(alphaRhoPhi, k_)
270  - fvm::laplacian(alpha*rho*DkEff(), k_)
271  ==
272  alpha()*rho()*G
273  - fvm::SuSp((2.0/3.0)*alpha()*rho()*divU, k_)
274  - fvm::Sp(Cmu_*alpha()*rho()*omega_(), k_)
275  + kSource()
276  + fvOptions(alpha, rho, k_)
277  );
278 
279  kEqn.ref().relax();
280  fvOptions.constrain(kEqn.ref());
281  solve(kEqn);
282  fvOptions.correct(k_);
283  bound(k_, this->kMin_);
284 
285  correctNut();
286 }
287 
288 
289 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
290 
291 } // End namespace RASModels
292 } // End namespace Foam
293 
294 // ************************************************************************* //
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
fvMatrix< scalar > fvScalarMatrix
Definition: fvMatricesFwd.H:42
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
void clear() const
If object pointer points to valid object:
Definition: tmpI.H:230
BasicTurbulenceModel::transportModel transportModel
Definition: kOmega.H:107
fv::options & fvOptions
surfaceScalarField & phi
kOmega(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.
Definition: kOmega.C:83
void correct(GeometricField< Type, fvPatchField, volMesh > &field)
Apply correction to field.
const dimensionedScalar G
Newtonian constant of gravitation.
T & ref() const
Return non-const reference or generate a fatal error.
Definition: tmpI.H:174
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
Generic dimensioned Type class.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
BasicTurbulenceModel::alphaField alphaField
Definition: kOmega.H:105
Eddy viscosity turbulence model base class.
Definition: eddyViscosity.H:49
Finite-volume options.
Definition: fvOptions.H:52
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
Templated abstract base class for RAS turbulence models.
Definition: RASModel.H:49
const dimensionSet dimVolume(pow3(dimLength))
Definition: dimensionSets.H:58
bool read(const char *, int32_t &)
Definition: int32IO.C:85
const Internal & v() const
Return a const-reference to the dimensioned internal field.
Definition: volFieldsI.H:31
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
BasicTurbulenceModel::rhoField rhoField
Definition: kOmega.H:106
A class for handling words, derived from string.
Definition: word.H:59
static word groupName(Name name, const word &group)
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: kOmega.C:210
void constrain(fvMatrix< Type > &eqn)
Apply constraints to equation.
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
Bound the given scalar field if it has gone unbounded.
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(rho0+psi *p, rhoMin);# 1 "/home/ubuntu/OpenFOAM-7/applications/solvers/multiphase/cavitatingFoam/alphavPsi.H" 1{ alphav=max(min((rho - rholSat)/(rhovSat - rholSat), scalar(1)), scalar(0));alphal=1.0 - alphav;Info<< "max-min alphav: "<< max(alphav).value()<< " "<< min(alphav).value()<< endl;psiModel-> correct()
Definition: pEqn.H:68
virtual void correctNut()
Definition: kOmega.C:40
virtual bool read()
Read RASProperties dictionary.
Definition: kOmega.C:190
tmp< surfaceScalarField > absolute(const tmp< surfaceScalarField > &tphi, const volVectorField &U)
Return the given relative flux in absolute form.
Definition: fvcMeshPhi.C:188
U
Definition: pEqn.H:72
volScalarField & bound(volScalarField &, const dimensionedScalar &lowerBound)
Bound the given scalar field if it has gone unbounded.
Definition: bound.C:33
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
fileType type(const fileName &, const bool checkVariants=true, const bool followLink=true)
Return the file type: directory or file.
Definition: POSIX.C:488
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
virtual tmp< fvScalarMatrix > kSource() const
Definition: kOmega.C:51
SolverPerformance< Type > solve(fvMatrix< Type > &, const word &)
Solve returning the solution statistics given convergence tolerance.
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
zeroField divU
Definition: alphaSuSp.H:3
A class for managing temporary objects.
Definition: PtrList.H:53
static options & New(const fvMesh &mesh)
Construct fvOptions and register to datbase if not present.
Definition: fvOptions.C:101
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:92
scalar nut
Namespace for OpenFOAM.
virtual tmp< fvScalarMatrix > omegaSource() const
Definition: kOmega.C:66