ReynoldsStress.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2015 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 "ReynoldsStress.H"
27 #include "fvc.H"
28 #include "fvm.H"
29 #include "wallFvPatch.H"
30 
31 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
32 
33 template<class BasicTurbulenceModel>
35 (
37 ) const
38 {
39  scalar kMin = this->kMin_.value();
40 
41  R.max
42  (
44  (
45  "zero",
46  R.dimensions(),
48  (
49  kMin, -GREAT, -GREAT,
50  kMin, -GREAT,
51  kMin
52  )
53  )
54  );
55 }
56 
57 
58 template<class BasicTurbulenceModel>
60 (
62 ) const
63 {
64  const fvPatchList& patches = this->mesh_.boundary();
65 
66  forAll(patches, patchi)
67  {
68  const fvPatch& curPatch = patches[patchi];
69 
70  if (isA<wallFvPatch>(curPatch))
71  {
73 
74  const scalarField& nutw = this->nut_.boundaryField()[patchi];
75 
76  const vectorField snGradU
77  (
78  this->U_.boundaryField()[patchi].snGrad()
79  );
80 
81  const vectorField& faceAreas
82  = this->mesh_.Sf().boundaryField()[patchi];
83 
84  const scalarField& magFaceAreas
85  = this->mesh_.magSf().boundaryField()[patchi];
86 
87  forAll(curPatch, facei)
88  {
89  // Calculate near-wall velocity gradient
90  tensor gradUw
91  = (faceAreas[facei]/magFaceAreas[facei])*snGradU[facei];
92 
93  // Calculate near-wall shear-stress tensor
94  tensor tauw = -nutw[facei]*2*dev(symm(gradUw));
95 
96  // Reset the shear components of the stress tensor
97  Rw[facei].xy() = tauw.xy();
98  Rw[facei].xz() = tauw.xz();
99  Rw[facei].yz() = tauw.yz();
100  }
101  }
102  }
103 }
104 
105 
106 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
107 
108 template<class BasicTurbulenceModel>
110 (
111  const word& modelName,
112  const alphaField& alpha,
113  const rhoField& rho,
114  const volVectorField& U,
115  const surfaceScalarField& alphaRhoPhi,
116  const surfaceScalarField& phi,
117  const transportModel& transport,
118  const word& propertiesName
119 )
120 :
121  BasicTurbulenceModel
122  (
123  modelName,
124  alpha,
125  rho,
126  U,
127  alphaRhoPhi,
128  phi,
129  transport,
130  propertiesName
131  ),
132 
133  couplingFactor_
134  (
136  (
137  "couplingFactor",
138  this->coeffDict_,
139  0.0
140  )
141  ),
142 
143  R_
144  (
145  IOobject
146  (
147  IOobject::groupName("R", U.group()),
148  this->runTime_.timeName(),
149  this->mesh_,
150  IOobject::MUST_READ,
151  IOobject::AUTO_WRITE
152  ),
153  this->mesh_
154  ),
155 
156  nut_
157  (
158  IOobject
159  (
160  IOobject::groupName("nut", U.group()),
161  this->runTime_.timeName(),
162  this->mesh_,
163  IOobject::MUST_READ,
164  IOobject::AUTO_WRITE
165  ),
166  this->mesh_
167  )
168 {
169  if (couplingFactor_.value() < 0.0 || couplingFactor_.value() > 1.0)
170  {
172  (
173  "ReynoldsStress::ReynoldsStress"
174  ) << "couplingFactor = " << couplingFactor_
175  << " is not in range 0 - 1" << nl
176  << exit(FatalError);
177  }
178 }
179 
180 
181 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
182 
183 template<class BasicTurbulenceModel>
185 {
187 }
188 
189 
190 template<class BasicTurbulenceModel>
193 {
194  return R_;
195 }
196 
197 
198 template<class BasicTurbulenceModel>
201 {
202  tmp<Foam::volScalarField> tk(0.5*tr(R_));
203  tk().rename("k");
204  return tk;
205 }
206 
207 
208 template<class BasicTurbulenceModel>
211 {
213  (
215  (
216  IOobject
217  (
218  IOobject::groupName("devRhoReff", this->U_.group()),
219  this->runTime_.timeName(),
220  this->mesh_,
221  IOobject::NO_READ,
222  IOobject::NO_WRITE
223  ),
224  this->alpha_*this->rho_*R_
225  - (this->alpha_*this->rho_*this->nu())
226  *dev(twoSymm(fvc::grad(this->U_)))
227  )
228  );
229 }
230 
231 
232 template<class BasicTurbulenceModel>
235 (
236  volVectorField& U
237 ) const
238 {
239  if (couplingFactor_.value() > 0.0)
240  {
241  return
242  (
243  fvc::div
244  (
245  this->alpha_*this->rho_*R_
246  + couplingFactor_
247  *this->alpha_*this->rho_*this->nut()*fvc::grad(U),
248  "div(devRhoReff)"
249  )
251  (
252  (1.0 - couplingFactor_)*this->alpha_*this->rho_*this->nut(),
253  U,
254  "laplacian(nuEff,U)"
255  )
256  - fvm::laplacian(this->alpha_*this->rho_*this->nuEff(), U)
257  - fvc::div(this->alpha_*this->rho_*this->nu()*dev2(T(fvc::grad(U))))
258  );
259  }
260  else
261  {
262  return
263  (
264  fvc::div(this->alpha_*this->rho_*R_)
266  (
267  this->alpha_*this->rho_*this->nut(),
268  U,
269  "laplacian(nuEff,U)"
270  )
271  - fvm::laplacian(this->alpha_*this->rho_*this->nuEff(), U)
272  - fvc::div(this->alpha_*this->rho_*this->nu()*dev2(T(fvc::grad(U))))
273  );
274  }
275 
276  return
277  (
278  - fvm::laplacian(this->alpha_*this->rho_*this->nuEff(), U)
279  - fvc::div((this->alpha_*this->rho_*this->nuEff())*dev2(T(fvc::grad(U))))
280  );
281 }
282 
283 
284 template<class BasicTurbulenceModel>
287 (
288  const volScalarField& rho,
289  volVectorField& U
290 ) const
291 {
292  return
293  (
294  - fvm::laplacian(this->alpha_*rho*this->nuEff(), U)
295  - fvc::div((this->alpha_*rho*this->nuEff())*dev2(T(fvc::grad(U))))
296  );
297 }
298 
299 
300 template<class BasicTurbulenceModel>
302 {
304 }
305 
306 
307 // ************************************************************************* //
void correctWallShearStress(volSymmTensorField &R) const
const Cmpt & xy() const
Definition: TensorI.H:167
virtual tmp< volScalarField > k() const
Return the turbulence kinetic energy.
virtual tmp< volSymmTensorField > devRhoReff() const
Return the effective stress tensor.
BasicTurbulenceModel::alphaField alphaField
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:61
GeometricBoundaryField & boundaryField()
Return reference to GeometricBoundaryField.
virtual void correct()=0
Solve the turbulence equations and correct the turbulence viscosity.
virtual tmp< fvVectorMatrix > divDevRhoReff(volVectorField &U) const
Return the source term for the momentum equation.
virtual bool read()=0
Re-read model coefficients if they have changed.
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/dm2/henry/OpenFOAM/OpenFOAM-3.0.x/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:74
A class for handling words, derived from string.
Definition: word.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
dimensionedScalar tr(const dimensionedSphericalTensor &dt)
SymmTensor< scalar > symmTensor
SymmTensor of scalars.
Definition: symmTensor.H:48
tmp< GeometricField< Type, fvPatchField, volMesh > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvcLaplacian.C:45
dimensionedSymmTensor dev2(const dimensionedSymmTensor &dt)
BasicTurbulenceModel::rhoField rhoField
const volScalarField & T
Definition: createFields.H:25
dimensionedSymmTensor dev(const dimensionedSymmTensor &dt)
patches[0]
BasicTurbulenceModel::transportModel transportModel
scalar tauw
Generic dimensioned Type class.
static const char nl
Definition: Ostream.H:260
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
scalar nut
#define forAll(list, i)
Definition: UList.H:421
label patchi
GeometricField< symmTensor, fvPatchField, volMesh > volSymmTensorField
Definition: volFieldsFwd.H:58
const dimensionSet & dimensions() const
Return dimensions.
void boundNormalStress(volSymmTensorField &R) const
virtual tmp< volSymmTensorField > R() const
Return the Reynolds stress tensor.
word group() const
Return group (extension part of name)
Definition: IOobject.C:263
#define FatalErrorIn(functionName)
Report an error message using Foam::FatalError.
Definition: error.H:314
dimensionedSymmTensor twoSymm(const dimensionedSymmTensor &dt)
error FatalError
void max(const dimensioned< Type > &)
bool read(const char *, int32_t &)
Definition: int32IO.C:87
volScalarField & nu
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:61
tmp< GeometricField< Type, fvPatchField, volMesh > > div(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcDiv.C:47
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
dimensionedSymmTensor symm(const dimensionedSymmTensor &dt)
const Cmpt & yz() const
Definition: TensorI.H:195
const Cmpt & xz() const
Definition: TensorI.H:174
A class for managing temporary objects.
Definition: PtrList.H:118
ReynoldsStress(const word &modelName, const alphaField &alpha, const rhoField &rho, const volVectorField &U, const surfaceScalarField &alphaRhoPhi, const surfaceScalarField &phi, const transportModel &transport, const word &propertiesName)
Construct from components.