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 | Copyright (C) 2011-2014 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 "GradientDispersionRAS.H"
27 #include "demandDrivenData.H"
28 #include "fvcGrad.H"
29 
30 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
31 
32 template<class CloudType>
34 (
35  const dictionary& dict,
36  CloudType& owner
37 )
38 :
39  DispersionRASModel<CloudType>(dict, owner),
40  gradkPtr_(NULL),
41  ownGradK_(false)
42 {}
43 
44 
45 template<class CloudType>
47 (
48  GradientDispersionRAS<CloudType>& dm
49 )
50 :
51  DispersionRASModel<CloudType>(dm),
52  gradkPtr_(dm.gradkPtr_),
53  ownGradK_(dm.ownGradK_)
54 {
55  dm.ownGradK_ = false;
56 }
57 
58 
59 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
60 
61 template<class CloudType>
63 {
64  cacheFields(false);
65 }
66 
67 
68 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
69 
70 template<class CloudType>
72 {
74 
75  if (store)
76  {
77  gradkPtr_ = fvc::grad(*this->kPtr_).ptr();
78  ownGradK_ = true;
79  }
80  else
81  {
82  if (ownGradK_)
83  {
84  deleteDemandDrivenData(gradkPtr_);
85  gradkPtr_ = NULL;
86  ownGradK_ = false;
87  }
88  }
89 }
90 
91 
92 template<class CloudType>
94 (
95  const scalar dt,
96  const label cellI,
97  const vector& U,
98  const vector& Uc,
99  vector& UTurb,
100  scalar& tTurb
101 )
102 {
103  cachedRandom& rnd = this->owner().rndGen();
104 
105  const scalar cps = 0.16432;
106 
107  const scalar k = this->kPtr_->internalField()[cellI];
108  const scalar epsilon =
109  this->epsilonPtr_->internalField()[cellI] + ROOTVSMALL;
110  const vector& gradk = this->gradkPtr_->internalField()[cellI];
111 
112  const scalar UrelMag = mag(U - Uc - UTurb);
113 
114  const scalar tTurbLoc =
115  min(k/epsilon, cps*pow(k, 1.5)/epsilon/(UrelMag + SMALL));
116 
117 
118  // Parcel is perturbed by the turbulence
119  if (dt < tTurbLoc)
120  {
121  tTurb += dt;
122 
123  if (tTurb > tTurbLoc)
124  {
125  tTurb = 0.0;
126 
127  scalar sigma = sqrt(2.0*k/3.0);
128  vector dir = -gradk/(mag(gradk) + SMALL);
129 
130  // Numerical Recipes... Ch. 7. Random Numbers...
131  scalar x1 = 0.0;
132  scalar x2 = 0.0;
133  scalar rsq = 10.0;
134  while ((rsq > 1.0) || (rsq == 0.0))
135  {
136  x1 = 2.0*rnd.sample01<scalar>() - 1.0;
137  x2 = 2.0*rnd.sample01<scalar>() - 1.0;
138  rsq = x1*x1 + x2*x2;
139  }
140 
141  scalar fac = sqrt(-2.0*log(rsq)/rsq);
142 
143  // In 2D calculations the -grad(k) is always
144  // away from the axis of symmetry
145  // This creates a 'hole' in the spray and to
146  // prevent this we let x1 be both negative/positive
147  if (this->owner().mesh().nSolutionD() == 2)
148  {
149  fac *= x1;
150  }
151  else
152  {
153  fac *= mag(x1);
154  }
155 
156  UTurb = sigma*fac*dir;
157  }
158  }
159  else
160  {
161  tTurb = GREAT;
162  UTurb = vector::zero;
163  }
164 
165  return Uc + UTurb;
166 }
167 
168 
169 // ************************************************************************* //
virtual void cacheFields(const bool store)
Cache carrier fields.
DsmcCloud< dsmcParcel > CloudType
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
dictionary dict
Vector< scalar > vector
A scalar version of the templated Vector.
Definition: vector.H:49
dynamicFvMesh & mesh
void deleteDemandDrivenData(DataPtr &dataPtr)
virtual vector update(const scalar dt, const label cellI, const vector &U, const vector &Uc, vector &UTurb, scalar &tTurb)
Update (disperse particles)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
dimensionedScalar log(const dimensionedScalar &ds)
GradientDispersionRAS(const dictionary &dict, CloudType &owner)
Construct from components.
virtual void cacheFields(const bool store)
Cache carrier fields.
Calculate the gradient of the given field.
Template functions to aid in the implementation of demand driven data.
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
tmp< GeometricField< typename outerProduct< vector, Type >::type, fvPatchField, volMesh >> grad(const GeometricField< Type, fvsPatchField, surfaceMesh > &ssf)
Definition: fvcGrad.C:52
virtual ~GradientDispersionRAS()
Destructor.