GradientDispersionRAS.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 | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2011-2016 OpenFOAM Foundation
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "GradientDispersionRAS.H"
29 #include "demandDrivenData.H"
30 #include "fvcGrad.H"
31 
32 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
34 template<class CloudType>
36 (
37  const dictionary& dict,
38  CloudType& owner
39 )
40 :
42  gradkPtr_(nullptr),
43  ownGradK_(false)
44 {}
45 
46 
47 template<class CloudType>
49 (
51 )
52 :
54  gradkPtr_(dm.gradkPtr_),
55  ownGradK_(dm.ownGradK_)
56 {
57  dm.ownGradK_ = false;
58 }
59 
60 
61 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
62 
63 template<class CloudType>
65 {
66  cacheFields(false);
67 }
68 
69 
70 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
71 
72 template<class CloudType>
74 {
76 
77  if (store)
78  {
79  gradkPtr_ = fvc::grad(*this->kPtr_).ptr();
80  ownGradK_ = true;
81  }
82  else
83  {
84  if (ownGradK_)
85  {
86  deleteDemandDrivenData(gradkPtr_);
87  gradkPtr_ = nullptr;
88  ownGradK_ = false;
89  }
90  }
91 }
92 
93 
94 template<class CloudType>
96 (
97  const scalar dt,
98  const label celli,
99  const vector& U,
100  const vector& Uc,
101  vector& UTurb,
102  scalar& tTurb
103 )
104 {
105  Random& rnd = this->owner().rndGen();
106 
107  const scalar cps = 0.16432;
108 
109  const scalar k = this->kPtr_->primitiveField()[celli];
110  const scalar epsilon =
111  this->epsilonPtr_->primitiveField()[celli] + ROOTVSMALL;
112  const vector& gradk = this->gradkPtr_->primitiveField()[celli];
113 
114  const scalar UrelMag = mag(U - Uc - UTurb);
115 
116  const scalar tTurbLoc =
117  min(k/epsilon, cps*pow(k, 1.5)/epsilon/(UrelMag + SMALL));
118 
119 
120  // Parcel is perturbed by the turbulence
121  if (dt < tTurbLoc)
122  {
123  tTurb += dt;
124 
125  if (tTurb > tTurbLoc)
126  {
127  tTurb = 0.0;
128 
129  const scalar sigma = sqrt(2.0*k/3.0);
130  const vector dir = -gradk/(mag(gradk) + SMALL);
131 
132  scalar fac = 0.0;
133 
134  // In 2D calculations the -grad(k) is always
135  // away from the axis of symmetry
136  // This creates a 'hole' in the spray and to
137  // prevent this we let fac be both negative/positive
138  if (this->owner().mesh().nSolutionD() == 2)
139  {
140  fac = rnd.GaussNormal<scalar>();
141  }
142  else
143  {
144  fac = mag(rnd.GaussNormal<scalar>());
145  }
146 
147  UTurb = sigma*fac*dir;
148  }
149  }
150  else
151  {
152  tTurb = GREAT;
153  UTurb = Zero;
154  }
155 
156  return Uc + UTurb;
157 }
158 
159 
160 // ************************************************************************* //
Foam::Random
Random number generator.
Definition: Random.H:59
Foam::GradientDispersionRAS::update
virtual vector update(const scalar dt, const label celli, const vector &U, const vector &Uc, vector &UTurb, scalar &tTurb)
Update (disperse particles)
Definition: GradientDispersionRAS.C:96
Foam::Random::GaussNormal
Type GaussNormal()
Definition: RandomTemplates.C:49
Foam::GradientDispersionRAS::~GradientDispersionRAS
virtual ~GradientDispersionRAS()
Destructor.
Definition: GradientDispersionRAS.C:64
Foam::GradientDispersionRAS::cacheFields
virtual void cacheFields(const bool store)
Cache carrier fields.
Definition: GradientDispersionRAS.C:73
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
demandDrivenData.H
Template functions to aid in the implementation of demand driven data.
Foam::GradientDispersionRAS::ownGradK_
bool ownGradK_
Take ownership of the grad(k)
Definition: GradientDispersionRAS.H:64
Foam::fac::grad
tmp< GeometricField< typename outerProduct< vector, Type >::type, faPatchField, areaMesh >> grad(const GeometricField< Type, faePatchField, edgeMesh > &ssf)
Definition: facGrad.C:56
Foam::min
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
Foam::DSMCCloud::rndGen
Random & rndGen()
Return reference to the random object.
Definition: DSMCCloudI.H:124
Foam::deleteDemandDrivenData
void deleteDemandDrivenData(DataPtr &dataPtr)
Definition: demandDrivenData.H:42
Foam::DSMCCloud::mesh
const fvMesh & mesh() const
Return reference to the mesh.
Definition: DSMCCloudI.H:44
Foam::pow
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Definition: dimensionedScalar.C:75
Foam::DSMCCloud
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:71
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
GradientDispersionRAS.H
Foam::GradientDispersionRAS::gradkPtr_
const volVectorField * gradkPtr_
Gradient of k.
Definition: GradientDispersionRAS.H:61
fac
Calculate the second temporal derivative.
Foam::polyMesh::nSolutionD
label nSolutionD() const
Return the number of valid solved-for dimensions in the mesh.
Definition: polyMesh.C:886
U
U
Definition: pEqn.H:72
Foam::Vector< scalar >
Foam::sqrt
dimensionedScalar sqrt(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:144
fvcGrad.H
Calculate the gradient of the given field.
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
k
label k
Boltzmann constant.
Definition: LISASMDCalcMethod2.H:41
epsilon
scalar epsilon
Definition: evaluateNearWall.H:7
Foam::GradientDispersionRAS
The velocity is perturbed in the direction of -grad(k), with a Gaussian random number distribution wi...
Definition: GradientDispersionRAS.H:50
sigma
dimensionedScalar sigma("sigma", dimMass/sqr(dimTime), transportProperties)
Foam::GradientDispersionRAS::GradientDispersionRAS
GradientDispersionRAS(const dictionary &dict, CloudType &owner)
Construct from components.
Definition: GradientDispersionRAS.C:36
Foam::DispersionRASModel
Base class for particle dispersion models based on RAS turbulence.
Definition: DispersionRASModel.H:49