GradientDispersionRAS.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) 2011-2024 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 #include "standardNormal.H"
30 
31 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
32 
33 template<class CloudType>
35 (
36  const dictionary& dict,
37  CloudType& owner
38 )
39 :
41  gradkPtr_(nullptr),
42  ownGradK_(false)
43 {}
44 
45 
46 template<class CloudType>
48 (
50 )
51 :
53  gradkPtr_(dm.gradkPtr_),
54  ownGradK_(dm.ownGradK_)
55 {
56  dm.ownGradK_ = false;
57 }
58 
59 
60 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
61 
62 template<class CloudType>
64 {
65  cacheFields(false);
66 }
67 
68 
69 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
70 
71 template<class CloudType>
73 {
75 
76  if (store)
77  {
78  gradkPtr_ = fvc::grad(*this->kPtr_).ptr();
79  ownGradK_ = true;
80  }
81  else
82  {
83  if (ownGradK_)
84  {
85  deleteDemandDrivenData(gradkPtr_);
86  gradkPtr_ = nullptr;
87  ownGradK_ = false;
88  }
89  }
90 }
91 
92 
93 template<class CloudType>
95 (
96  const scalar dt,
97  const label celli,
98  const vector& U,
99  const vector& Uc,
100  vector& UTurb,
101  scalar& tTurb
102 )
103 {
104  distributions::standardNormal& stdNormal = this->owner().stdNormal();
105 
106  const scalar cps = 0.16432;
107 
108  const scalar k = this->kPtr_->primitiveField()[celli];
109  const scalar epsilon =
110  this->epsilonPtr_->primitiveField()[celli] + rootVSmall;
111  const vector& gradk = this->gradkPtr_->primitiveField()[celli];
112 
113  const scalar UrelMag = mag(U - Uc - UTurb);
114 
115  const scalar tTurbLoc =
116  min(k/epsilon, cps*pow(k, 1.5)/epsilon/(UrelMag + small));
117 
118 
119  // Parcel is perturbed by the turbulence
120  if (dt < tTurbLoc)
121  {
122  tTurb += dt;
123 
124  if (tTurb > tTurbLoc)
125  {
126  tTurb = 0.0;
127 
128  const scalar sigma = sqrt(2.0*k/3.0);
129  const vector dir = -gradk/(mag(gradk) + small);
130 
131  scalar fac = 0.0;
132 
133  // In 2D calculations the -grad(k) is always
134  // away from the axis of symmetry
135  // This creates a 'hole' in the spray and to
136  // prevent this we let fac be both negative/positive
137  if (this->owner().mesh().nSolutionD() == 2)
138  {
139  fac = stdNormal.sample();
140  }
141  else
142  {
143  fac = mag(stdNormal.sample());
144  }
145 
146  UTurb = sigma*fac*dir;
147  }
148  }
149  else
150  {
151  tTurb = great;
152  UTurb = Zero;
153  }
154 
155  return Uc + UTurb;
156 }
157 
158 
159 // ************************************************************************* //
label k
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:80
Base class for particle dispersion models based on RAS turbulence.
virtual void cacheFields(const bool store)
Cache carrier fields.
The velocity is perturbed in the direction of -grad(k), with a Gaussian random number distribution wi...
virtual ~GradientDispersionRAS()
Destructor.
virtual void cacheFields(const bool store)
Cache carrier fields.
GradientDispersionRAS(const dictionary &dict, CloudType &owner)
Construct from components.
virtual vector update(const scalar dt, const label celli, const vector &U, const vector &Uc, vector &UTurb, scalar &tTurb)
Update (disperse particles)
bool ownGradK_
Take ownership of the grad(k)
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
Standard normal distribution. Not selectable.
virtual scalar sample() const
Sample the distribution.
Template functions to aid in the implementation of demand driven data.
const scalar epsilon
Calculate the gradient of the given field.
U
Definition: pEqn.H:72
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m^2/K^4].
tmp< VolField< typename outerProduct< vector, Type >::type > > grad(const SurfaceField< Type > &ssf)
Definition: fvcGrad.C:46
static const zero Zero
Definition: zero.H:97
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
void deleteDemandDrivenData(DataType *&dataPtr)
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensionedScalar sqrt(const dimensionedScalar &ds)
dimensioned< scalar > mag(const dimensioned< Type > &)
dictionary dict