FreeFOAM The Cross-Platform CFD Toolkit
VariableHardSphere.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) 2009-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 "VariableHardSphere.H"
27 
28 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
29 
30 template <class CloudType>
32 (
33  const dictionary& dict,
34  CloudType& cloud
35 )
36 :
37  BinaryCollisionModel<CloudType>(dict, cloud, typeName),
38  Tref_(readScalar(this->coeffDict().lookup("Tref")))
39 {}
40 
41 
42 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
43 
44 template <class CloudType>
46 {}
47 
48 
49 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
50 
51 
52 template <class CloudType>
54 (
55  label typeIdP,
56  label typeIdQ,
57  const vector& UP,
58  const vector& UQ
59 ) const
60 {
61  const CloudType& cloud(this->owner());
62 
63  scalar dPQ =
64  0.5
65  *(
66  cloud.constProps(typeIdP).d()
67  + cloud.constProps(typeIdQ).d()
68  );
69 
70  scalar omegaPQ =
71  0.5
72  *(
73  cloud.constProps(typeIdP).omega()
74  + cloud.constProps(typeIdQ).omega()
75  );
76 
77  scalar cR = mag(UP - UQ);
78 
79  if (cR < VSMALL)
80  {
81  return 0;
82  }
83 
84  scalar mP = cloud.constProps(typeIdP).mass();
85 
86  scalar mQ = cloud.constProps(typeIdQ).mass();
87 
88  scalar mR = mP*mQ/(mP + mQ);
89 
90  // calculating cross section = pi*dPQ^2, where dPQ is from Bird, eq. 4.79
91  scalar sigmaTPQ =
93  *pow(2.0*CloudType::kb*Tref_/(mR*cR*cR), omegaPQ - 0.5)
94  /exp(Foam::lgamma(2.5 - omegaPQ));
95 
96  return sigmaTPQ*cR;
97 }
98 
99 
100 template <class CloudType>
102 (
103  label typeIdP,
104  label typeIdQ,
105  vector& UP,
106  vector& UQ,
107  scalar& EiP,
108  scalar& EiQ
109 )
110 {
111  CloudType& cloud(this->owner());
112 
113  Random& rndGen(cloud.rndGen());
114 
115  scalar mP = cloud.constProps(typeIdP).mass();
116 
117  scalar mQ = cloud.constProps(typeIdQ).mass();
118 
119  vector Ucm = (mP*UP + mQ*UQ)/(mP + mQ);
120 
121  scalar cR = mag(UP - UQ);
122 
123  scalar cosTheta = 2.0*rndGen.scalar01() - 1.0;
124 
125  scalar sinTheta = sqrt(1.0 - cosTheta*cosTheta);
126 
127  scalar phi = 2.0*mathematicalConstant::pi*rndGen.scalar01();
128 
129  vector postCollisionRelU =
130  cR
131  *vector
132  (
133  cosTheta,
134  sinTheta*cos(phi),
135  sinTheta*sin(phi)
136  );
137 
138  UP = Ucm + postCollisionRelU*mQ/(mP + mQ);
139 
140  UQ = Ucm - postCollisionRelU*mP/(mP + mQ);
141 }
142 
143 
144 // ************************ vim: set sw=4 sts=4 et: ************************ //