PDRkEpsilon.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-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 "PDRkEpsilon.H"
27 #include "PDRDragModel.H"
29 
30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
34 namespace compressible
35 {
36 namespace RASModels
37 {
38 
39 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
40 
41 defineTypeNameAndDebug(PDRkEpsilon, 0);
42 addToRunTimeSelectionTable(RASModel, PDRkEpsilon, dictionary);
43 
44 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
45 
47 (
48  const geometricOneField& alpha,
49  const volScalarField& rho,
50  const volVectorField& U,
51  const surfaceScalarField& alphaRhoPhi,
52  const surfaceScalarField& phi,
53  const fluidThermo& thermophysicalModel,
54  const word& turbulenceModelName,
55  const word& modelName
56 )
57 :
58  Foam::RASModels::kEpsilon<EddyDiffusivity<compressible::turbulenceModel>>
59  (
60  geometricOneField(),
61  rho,
62  U,
63  phi,
64  phi,
65  thermophysicalModel,
66  turbulenceModelName,
67  modelName
68  ),
69 
70  C4_
71  (
72  dimensioned<scalar>::lookupOrAddToDict
73  (
74  "C4",
75  coeffDict_,
76  0.1
77  )
78  )
79 {}
80 
81 
82 // * * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * //
83 
84 PDRkEpsilon::~PDRkEpsilon()
85 {}
86 
87 
88 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
89 
90 bool PDRkEpsilon::read()
91 {
92  if (RASModel::read())
93  {
94  C4_.readIfPresent(coeffDict_);
95  return true;
96  }
97  else
98  {
99  return false;
100  }
101 }
102 
103 
105 {
106  if (!turbulence_)
107  {
108  // Re-calculate viscosity
109  nut_ = Cmu_*sqr(k_)/epsilon_;
110  nut_.correctBoundaryConditions();
111 
112  // Re-calculate thermal diffusivity
113  //***HGWalphat_ = mut_/Prt_;
114  // alphat_.correctBoundaryConditions();
115 
116  return;
117  }
118 
120 
122 
123  if (mesh_.moving())
124  {
125  divU += fvc::div(mesh_.phi());
126  }
127 
128  tmp<volTensorField> tgradU = fvc::grad(U_);
129  volScalarField G(GName(), rho_*nut_*(tgradU() && dev(twoSymm(tgradU()))));
130  tgradU.clear();
131 
132  // Update espsilon and G at the wall
133  epsilon_.boundaryFieldRef().updateCoeffs();
134 
135  // Add the blockage generation term so that it is included consistently
136  // in both the k and epsilon equations
137  const volScalarField& betav =
138  U_.db().lookupObject<volScalarField>("betav");
139 
140  const volScalarField& Lobs =
141  U_.db().lookupObject<volScalarField>("Lobs");
142 
143  const PDRDragModel& drag =
144  U_.db().lookupObject<PDRDragModel>("PDRDragModel");
145 
146  volScalarField GR(drag.Gk());
147 
148  volScalarField LI
149  (C4_*(Lobs + dimensionedScalar("minLength", dimLength, rootVSmall)));
150 
151  // Dissipation equation
152  tmp<fvScalarMatrix> epsEqn
153  (
154  betav*fvm::ddt(rho_, epsilon_)
155  + fvm::div(phi_, epsilon_)
156  - fvm::laplacian(rho_*DepsilonEff(), epsilon_)
157  ==
158  C1_*betav*G*epsilon_/k_
159  + 1.5*pow(Cmu_, 3.0/4.0)*GR*sqrt(k_)/LI
160  - fvm::SuSp(((2.0/3.0)*C1_)*betav*rho_*divU, epsilon_)
161  - fvm::Sp(C2_*betav*rho_*epsilon_/k_, epsilon_)
162  );
163 
164  epsEqn.ref().relax();
165 
166  epsEqn.ref().boundaryManipulate(epsilon_.boundaryFieldRef());
167 
168  solve(epsEqn);
169  bound(epsilon_, epsilonMin_);
170 
171 
172  // Turbulent kinetic energy equation
173 
174  tmp<fvScalarMatrix> kEqn
175  (
176  betav*fvm::ddt(rho_, k_)
177  + fvm::div(phi_, k_)
178  - fvm::laplacian(rho_*DkEff(), k_)
179  ==
180  betav*G + GR
181  - fvm::SuSp((2.0/3.0)*betav*rho_*divU, k_)
182  - fvm::Sp(betav*rho_*epsilon_/k_, k_)
183  );
184 
185  kEqn.ref().relax();
186  solve(kEqn);
187  bound(k_, kMin_);
188 
189  // Re-calculate viscosity
190  nut_ = Cmu_*sqr(k_)/epsilon_;
191  nut_.correctBoundaryConditions();
192 
193  // Re-calculate thermal diffusivity
194  //***HGWalphat_ = mut_/Prt_;
195  // alphat_.correctBoundaryConditions();
196 }
197 
198 
199 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
200 
201 } // End namespace RASModels
202 } // End namespace compressible
203 } // End namespace Foam
204 
205 // ************************************************************************* //
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
tmp< fvMatrix< Type > > SuSp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
surfaceScalarField & phi
const dimensionedScalar G
Newtonian constant of gravitation.
dimensionedSymmTensor sqr(const dimensionedVector &dv)
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
PDRkEpsilon(const geometricOneField &alpha, const volScalarField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const fluidThermo &thermophysicalModel, const word &turbulenceModelName=turbulenceModel::typeName, const word &modelName=typeName)
Construct from components.
ThermalDiffusivity< CompressibleTurbulenceModel< fluidThermo > > turbulenceModel
dimensionedScalar sqrt(const dimensionedScalar &ds)
void read(Istream &, label &, const dictionary &)
In-place read with dictionary lookup.
tmp< fvMatrix< Type > > Sp(const volScalarField::Internal &, const GeometricField< Type, fvPatchField, volMesh > &)
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
virtual bool read()
Read model coefficients if they have changed.
Definition: RASModel.C:169
const Type & lookupObject(const word &name) const
Lookup and return the object of the given Type.
const compressible::RASModel & turbulence_
Definition: PDRDragModel.H:64
Macros for easy insertion into run-time selection tables.
const volScalarField & rho_
Definition: PDRDragModel.H:65
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
const volVectorField & U_
Definition: PDRDragModel.H:66
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
SolverPerformance< Type > solve(fvMatrix< Type > &, const dictionary &)
Solve returning the solution statistics given convergence tolerance.
tmp< fvMatrix< Type > > ddt(const GeometricField< Type, fvPatchField, volMesh > &vf)
Definition: fvmDdt.C:46
const surfaceScalarField & phi_
Definition: PDRDragModel.H:67
addToRunTimeSelectionTable(ensightPart, ensightPartCells, istream)
tmp< fvMatrix< Type > > div(const surfaceScalarField &flux, const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmDiv.C:46
Info<< "Predicted p max-min : "<< max(p).value()<< " "<< min(p).value()<< endl;rho==max(psi *p+alphal *rhol0+((alphav *psiv+alphal *psil) - psi) *pSat, rhoMin);# 1 "/home/ubuntu/OpenFOAM-6/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:72
defineTypeNameAndDebug(combustionModel, 0)
RASModel< EddyDiffusivity< turbulenceModel > > RASModel
static tmp< GeometricField< Type, fvsPatchField, surfaceMesh > > interpolate(const GeometricField< Type, fvPatchField, volMesh > &tvf, const surfaceScalarField &faceFlux, Istream &schemeData)
Interpolate field onto faces using scheme given by Istream.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
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
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
zeroField divU
Definition: alphaSuSp.H:3
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:361
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
virtual void correct()
Solve the turbulence equations and correct the turbulence viscosity.
Definition: RASModel.C:192
Namespace for OpenFOAM.