FreeFOAM The Cross-Platform CFD Toolkit
stochasticDispersionRAS.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) 1991-2010 OpenCFD Ltd.
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 <OpenFOAM/error.H>
27 
30 
31 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
32 
33 namespace Foam
34 {
35 
36 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
37 
38 defineTypeNameAndDebug(stochasticDispersionRAS, 0);
39 
41 (
42  dispersionModel,
43  stochasticDispersionRAS,
44  dictionary
45 );
46 
47 
48 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
49 
50 // Construct from components
52 (
53  const dictionary& dict,
54  spray& sm
55 )
56 :
57  dispersionRASModel(dict, sm)
58 {}
59 
60 
61 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
62 
64 {}
65 
66 
67 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
68 
70 {
71 
72  const scalar cps = 0.16432;
73  const vector one(1.0, 1.0, 1.0);
74 
75  scalar dt = spray_.runTime().deltaT().value();
76  const volScalarField& k = turbulence().k();
77  //volVectorField gradk = fvc::grad(k);
78  const volScalarField& epsilon = turbulence().epsilon();
79  const volVectorField& U = spray_.U();
80 
81  for
82  (
84  elmnt != spray_.end();
85  ++elmnt
86  )
87  {
88  label celli = elmnt().cell();
89  scalar UrelMag = mag(elmnt().U() - U[celli] - elmnt().Uturb());
90 
91  scalar Tturb = min
92  (
93  k[celli]/epsilon[celli],
94  cps*pow(k[celli], 1.5)/epsilon[celli]/(UrelMag + SMALL)
95  );
96 
97  // parcel is perturbed by the turbulence
98  if (dt < Tturb)
99  {
100  elmnt().tTurb() += dt;
101 
102  if (elmnt().tTurb() > Tturb)
103  {
104  elmnt().tTurb() = 0.0;
105 
106  scalar sigma = sqrt(2.0*k[celli]/3.0);
107  vector dir = 2.0*spray_.rndGen().vector01() - one;
108  dir /= mag(dir) + SMALL;
109 
110  // numerical recipes... Ch. 7. Random Numbers...
111  scalar x1,x2;
112  scalar rsq = 10.0;
113  while((rsq > 1.0) || (rsq == 0.0))
114  {
115  x1 = 2.0*spray_.rndGen().scalar01() - 1.0;
116  x2 = 2.0*spray_.rndGen().scalar01() - 1.0;
117  rsq = x1*x1 + x2*x2;
118  }
119 
120  scalar fac = sqrt(-2.0*log(rsq)/rsq);
121 
122  fac *= mag(x1);
123 
124  elmnt().Uturb() = sigma*fac*dir;
125 
126  }
127  }
128  else
129  {
130  elmnt().tTurb() = GREAT;
131  elmnt().Uturb() = vector::zero;
132  }
133  }
134 }
135 
136 
137 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
138 
139 } // End namespace Foam
140 
141 // ************************ vim: set sw=4 sts=4 et: ************************ //