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  <
46  >
47  (
48  type,
49  alpha,
50  rho,
51  U,
52  alphaRhoPhi,
53  phi,
54  phase,
55  propertiesName
56  ),
57 
58  phase_(phase),
59 
60  alphaMax_(readScalar(coeffDict_.lookup("alphaMax"))),
61  preAlphaExp_(readScalar(coeffDict_.lookup("preAlphaExp"))),
62  expMax_(readScalar(coeffDict_.lookup("expMax"))),
63  g0_
64  (
65  "g0",
66  dimensionSet(1, -1, -2, 0, 0),
67  coeffDict_.lookup("g0")
68  )
69 {
70  nut_ == dimensionedScalar("zero", nut_.dimensions(), 0.0);
71 
72  if (type == typeName)
73  {
74  printCoeffs(type);
75  }
76 }
77 
78 
79 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
80 
82 {}
83 
84 
85 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
86 
88 {
89  if
90  (
91  eddyViscosity
92  <
93  RASModel<EddyDiffusivity<phaseCompressibleTurbulenceModel>>
94  >::read()
95  )
96  {
97  coeffDict().lookup("alphaMax") >> alphaMax_;
98  coeffDict().lookup("preAlphaExp") >> preAlphaExp_;
99  coeffDict().lookup("expMax") >> expMax_;
100  g0_.readIfPresent(coeffDict());
101 
102  return true;
103  }
104  else
105  {
106  return false;
107  }
108 }
109 
110 
113 {
115  return nut_;
116 }
117 
118 
121 {
123  return nut_;
124 }
125 
126 
129 {
130  return tmp<volSymmTensorField>
131  (
133  (
134  IOobject
135  (
136  IOobject::groupName("R", U_.group()),
137  runTime_.timeName(),
138  mesh_,
141  ),
142  mesh_,
143  dimensioned<symmTensor>
144  (
145  "R",
146  dimensionSet(0, 2, -2, 0, 0),
147  Zero
148  )
149  )
150  );
151 }
152 
153 
156 {
157  tmp<volScalarField> tpPrime
158  (
159  g0_
160  *min
161  (
162  exp(preAlphaExp_*(alpha_ - alphaMax_)),
163  expMax_
164  )
165  );
166 
167  volScalarField::Boundary& bpPrime =
168  tpPrime.ref().boundaryFieldRef();
169 
170  forAll(bpPrime, patchi)
171  {
172  if (!bpPrime[patchi].coupled())
173  {
174  bpPrime[patchi] == 0;
175  }
176  }
177 
178  return tpPrime;
179 }
180 
181 
184 {
185  tmp<surfaceScalarField> tpPrime
186  (
187  g0_
188  *min
189  (
190  exp(preAlphaExp_*(fvc::interpolate(alpha_) - alphaMax_)),
191  expMax_
192  )
193  );
194 
195  surfaceScalarField::Boundary& bpPrime =
196  tpPrime.ref().boundaryFieldRef();
197 
198  forAll(bpPrime, patchi)
199  {
200  if (!bpPrime[patchi].coupled())
201  {
202  bpPrime[patchi] == 0;
203  }
204  }
205 
206  return tpPrime;
207 }
208 
209 
212 {
213  return tmp<volSymmTensorField>
214  (
216  (
217  IOobject
218  (
219  IOobject::groupName("devRhoReff", U_.group()),
220  runTime_.timeName(),
221  mesh_,
224  ),
225  mesh_,
226  dimensioned<symmTensor>
227  (
228  "R",
229  rho_.dimensions()*dimensionSet(0, 2, -2, 0, 0),
230  Zero
231  )
232  )
233  );
234 }
235 
236 
239 (
240  volVectorField& U
241 ) const
242 {
243  return tmp<fvVectorMatrix>
244  (
245  new fvVectorMatrix
246  (
247  U,
248  rho_.dimensions()*dimensionSet(0, 4, -2, 0, 0)
249  )
250  );
251 }
252 
253 
255 {}
256 
257 
258 // ************************************************************************* //
#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
ThermalDiffusivity< PhaseCompressibleTurbulenceModel< phaseModel > > phaseCompressibleTurbulenceModel
Typedef for phaseCompressibleTurbulenceModel.
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;.