phasePressureModel.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 "phasePressureModel.H"
27 #include "twoPhaseSystem.H"
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 Foam::RASModels::phasePressureModel::phasePressureModel
32 (
33  const volScalarField& alpha,
34  const volScalarField& rho,
35  const volVectorField& U,
36  const surfaceScalarField& alphaRhoPhi,
37  const surfaceScalarField& phi,
38  const transportModel& phase,
39  const word& propertiesName,
40  const word& type
41 )
42 :
43  eddyViscosity
44  <
45  RASModel<EddyDiffusivity<ThermalDiffusivity
46  <
47  PhaseCompressibleTurbulenceModel<phaseModel>
48  >>>
49  >
50  (
51  type,
52  alpha,
53  rho,
54  U,
55  alphaRhoPhi,
56  phi,
57  phase,
58  propertiesName
59  ),
60 
61  phase_(phase),
62 
63  alphaMax_(readScalar(coeffDict_.lookup("alphaMax"))),
64  preAlphaExp_(readScalar(coeffDict_.lookup("preAlphaExp"))),
65  expMax_(readScalar(coeffDict_.lookup("expMax"))),
66  g0_
67  (
68  "g0",
69  dimensionSet(1, -1, -2, 0, 0),
70  coeffDict_.lookup("g0")
71  )
72 {
73  nut_ == dimensionedScalar("zero", nut_.dimensions(), 0.0);
74 
75  if (type == typeName)
76  {
77  printCoeffs(type);
78  }
79 }
80 
81 
82 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
83 
85 {}
86 
87 
88 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
89 
91 {
92  if
93  (
94  eddyViscosity
95  <
96  RASModel<EddyDiffusivity<ThermalDiffusivity
97  <
98  PhaseCompressibleTurbulenceModel<phaseModel>
99  >>>
100  >::read()
101  )
102  {
103  coeffDict().lookup("alphaMax") >> alphaMax_;
104  coeffDict().lookup("preAlphaExp") >> preAlphaExp_;
105  coeffDict().lookup("expMax") >> expMax_;
106  g0_.readIfPresent(coeffDict());
107 
108  return true;
109  }
110  else
111  {
112  return false;
113  }
114 }
115 
116 
119 {
121  return nut_;
122 }
123 
124 
127 {
129  return nut_;
130 }
131 
132 
135 {
136  return tmp<volSymmTensorField>
137  (
139  (
140  IOobject
141  (
142  IOobject::groupName("R", U_.group()),
143  runTime_.timeName(),
144  mesh_,
147  ),
148  mesh_,
149  dimensioned<symmTensor>
150  (
151  "R",
152  dimensionSet(0, 2, -2, 0, 0),
153  Zero
154  )
155  )
156  );
157 }
158 
159 
162 {
163  tmp<volScalarField> tpPrime
164  (
165  g0_
166  *min
167  (
168  exp(preAlphaExp_*(alpha_ - alphaMax_)),
169  expMax_
170  )
171  );
172 
173  volScalarField::Boundary& bpPrime =
174  tpPrime.ref().boundaryFieldRef();
175 
176  forAll(bpPrime, patchi)
177  {
178  if (!bpPrime[patchi].coupled())
179  {
180  bpPrime[patchi] == 0;
181  }
182  }
183 
184  return tpPrime;
185 }
186 
187 
190 {
191  tmp<surfaceScalarField> tpPrime
192  (
193  g0_
194  *min
195  (
196  exp(preAlphaExp_*(fvc::interpolate(alpha_) - alphaMax_)),
197  expMax_
198  )
199  );
200 
201  surfaceScalarField::Boundary& bpPrime =
202  tpPrime.ref().boundaryFieldRef();
203 
204  forAll(bpPrime, patchi)
205  {
206  if (!bpPrime[patchi].coupled())
207  {
208  bpPrime[patchi] == 0;
209  }
210  }
211 
212  return tpPrime;
213 }
214 
215 
218 {
219  return tmp<volSymmTensorField>
220  (
222  (
223  IOobject
224  (
225  IOobject::groupName("devRhoReff", U_.group()),
226  runTime_.timeName(),
227  mesh_,
230  ),
231  mesh_,
232  dimensioned<symmTensor>
233  (
234  "R",
235  rho_.dimensions()*dimensionSet(0, 2, -2, 0, 0),
236  Zero
237  )
238  )
239  );
240 }
241 
242 
245 (
246  volVectorField& U
247 ) const
248 {
249  return tmp<fvVectorMatrix>
250  (
251  new fvVectorMatrix
252  (
253  U,
254  rho_.dimensions()*dimensionSet(0, 4, -2, 0, 0)
255  )
256  );
257 }
258 
259 
261 {}
262 
263 
264 // ************************************************************************* //
#define readScalar
Definition: doubleScalar.C:38
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
Definition: volFieldsFwd.H:58
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
GeometricField< vector, fvPatchField, volMesh > volVectorField
Definition: volFieldsFwd.H:55
virtual tmp< volScalarField > pPrime() const
Return the phase-pressure&#39;.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
stressControl lookup("compactNormalStress") >> compactNormalStress
dimensionedScalar exp(const dimensionedScalar &ds)
static word groupName(Name name, const word &group)
virtual tmp< volSymmTensorField > R() const
Return the Reynolds stress tensor.
static const zero Zero
Definition: zero.H:97
virtual void correct()
Solve the kinetic theory equations and correct the viscosity.
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
virtual bool read()
Re-read model coefficients if they have changed.
RASModel< EddyDiffusivity< turbulenceModel > > RASModel
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
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.
virtual tmp< fvVectorMatrix > divDevRhoReff(volVectorField &U) const
Return the source term for the momentum equation.
virtual ~phasePressureModel()
Destructor.
label patchi
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
fvMatrix< vector > fvVectorMatrix
Definition: fvMatricesFwd.H:45
A class for managing temporary objects.
Definition: PtrList.H:53
GeometricField< scalar, fvsPatchField, surfaceMesh > surfaceScalarField
virtual tmp< volSymmTensorField > devRhoReff() const
Return the effective stress tensor.
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Definition: error.H:366
virtual tmp< volScalarField > epsilon() const
Return the turbulence kinetic energy dissipation rate.
virtual tmp< surfaceScalarField > pPrimef() const
Return the face-phase-pressure&#39;.