1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 1991-2008 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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 the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
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
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
25 \*---------------------------------------------------------------------------*/
29 #include "stochasticDispersionRAS.H"
30 #include "addToRunTimeSelectionTable.H"
32 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
37 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
39 defineTypeNameAndDebug(stochasticDispersionRAS, 0);
41 addToRunTimeSelectionTable
44 stochasticDispersionRAS,
49 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
51 // Construct from components
52 stochasticDispersionRAS::stochasticDispersionRAS
54 const dictionary& dict,
58 dispersionRASModel(dict, sm)
62 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
64 stochasticDispersionRAS::~stochasticDispersionRAS()
68 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
70 void stochasticDispersionRAS::disperseParcels() const
73 const scalar cps = 0.16432;
74 const vector one(1.0, 1.0, 1.0);
76 scalar dt = spray_.runTime().deltaT().value();
77 const volScalarField& k = turbulence().k();
78 //volVectorField gradk = fvc::grad(k);
79 const volScalarField& epsilon = turbulence().epsilon();
80 const volVectorField& U = spray_.U();
84 spray::iterator elmnt = spray_.begin();
85 elmnt != spray_.end();
89 label celli = elmnt().cell();
90 scalar UrelMag = mag(elmnt().U() - U[celli] - elmnt().Uturb());
94 k[celli]/epsilon[celli],
95 cps*pow(k[celli], 1.5)/epsilon[celli]/(UrelMag + SMALL)
98 // parcel is perturbed by the turbulence
101 elmnt().tTurb() += dt;
103 if (elmnt().tTurb() > Tturb)
105 elmnt().tTurb() = 0.0;
107 scalar sigma = sqrt(2.0*k[celli]/3.0);
108 vector dir = 2.0*spray_.rndGen().vector01() - one;
109 dir /= mag(dir) + SMALL;
111 // numerical recipes... Ch. 7. Random Numbers...
114 while((rsq > 1.0) || (rsq == 0.0))
116 x1 = 2.0*spray_.rndGen().scalar01() - 1.0;
117 x2 = 2.0*spray_.rndGen().scalar01() - 1.0;
121 scalar fac = sqrt(-2.0*log(rsq)/rsq);
125 elmnt().Uturb() = sigma*fac*dir;
131 elmnt().tTurb() = GREAT;
132 elmnt().Uturb() = vector::zero;
138 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
140 } // End namespace Foam
142 // ************************************************************************* //