TrajectoryCollision.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-2022 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 "TrajectoryCollision.H"
27 
28 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
29 
30 template<class CloudType>
32 (
33  typename CloudType::parcelType::trackingData& td
34 )
35 {
37 }
38 
39 
40 template<class CloudType>
42 (
43  const scalar dt,
44  parcelType& p1,
45  parcelType& p2,
46  scalar& m1,
47  scalar& m2
48 )
49 {
50  bool coalescence = false;
51 
52  const point pos1 = p1.position(this->owner().mesh());
53  const point pos2 = p2.position(this->owner().mesh());
54 
55  const vector& U1 = p1.U();
56  const vector& U2 = p2.U();
57 
58  vector URel = U1 - U2;
59 
60  vector d = pos2 - pos1;
61  scalar magd = mag(d);
62 
63  scalar vAlign = URel & (d/(magd + rootVSmall));
64 
65  if (vAlign > 0)
66  {
67  const scalar d1 = p1.d();
68  const scalar d2 = p2.d();
69 
70  scalar sumD = d1 + d2;
71 
72  if (vAlign*dt > magd - 0.5*sumD)
73  {
74  scalar magU1 = mag(U1) + rootVSmall;
75  scalar magU2 = mag(U2) + rootVSmall;
76  vector n1 = U1/magU1;
77  vector n2 = U2/magU2;
78 
79  scalar n1n2 = n1 & n2;
80  scalar n1d = n1 & d;
81  scalar n2d = n2 & d;
82 
83  scalar det = 1.0 - sqr(n1n2);
84 
85  scalar alpha = great;
86  scalar beta = great;
87 
88  if (mag(det) > 1.0e-4)
89  {
90  beta = -(n2d - n1n2*n1d)/det;
91  alpha = n1d + n1n2*beta;
92  }
93 
94  alpha /= magU1*dt;
95  beta /= magU2*dt;
96 
97  // is collision possible within this timestep
98  if ((alpha > 0) && (alpha < 1.0) && (beta > 0) && (beta < 1.0))
99  {
100  vector p1c = pos1 + alpha*U1*dt;
101  vector p2c = pos2 + beta*U2*dt;
102 
103  scalar closestDist = mag(p1c - p2c);
104 
105  scalar collProb =
106  pow(0.5*sumD/max(0.5*sumD, closestDist), cSpace_)
107  *exp(-cTime_*mag(alpha - beta));
108 
109  scalar xx = this->owner().rndGen().template sample01<scalar>();
110 
111  // collision occurs
112  if (xx > collProb)
113  {
114  if (d1 > d2)
115  {
116  coalescence = this->collideSorted(dt, p1, p2, m1, m2);
117  }
118  else
119  {
120  coalescence = this->collideSorted(dt, p2, p1, m2, m1);
121  }
122  }
123  }
124  }
125  }
126 
127  return coalescence;
128 }
129 
130 
131 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
132 
133 template<class CloudType>
135 (
136  const dictionary& dict,
137  CloudType& owner
138 )
139 :
140  ORourkeCollision<CloudType>(dict, owner, typeName),
141  cSpace_(this->coeffDict().template lookup<scalar>("cSpace")),
142  cTime_(this->coeffDict().template lookup<scalar>("cTime"))
143 {}
144 
145 
146 template<class CloudType>
148 (
150 )
151 :
153  cSpace_(cm.cSpace_),
154  cTime_(cm.cTime_)
155 {}
156 
157 
158 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
159 
160 template<class CloudType>
162 {}
163 
164 
165 // ************************************************************************* //
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:80
Collision model by P.J. O'Rourke.
CloudType::parcelType parcelType
Convenience typedef to the cloud's parcel type.
Trajectory collision model by N. Nordin, based on O'Rourke's collision model.
virtual void collide(typename CloudType::parcelType::trackingData &td)
Main collision routine.
virtual ~TrajectoryCollision()
Destructor.
virtual bool collideParcels(const scalar dt, parcelType &p1, parcelType &p2, scalar &m1, scalar &m2)
Collide parcels and return true if mass has changed.
TrajectoryCollision(const dictionary &dict, CloudType &cloud)
Construct from dictionary.
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
volScalarField alpha(IOobject("alpha", runTime.name(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
dimensionedScalar det(const dimensionedSphericalTensor &dt)
dimensionedScalar exp(const dimensionedScalar &ds)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensioned< scalar > mag(const dimensioned< Type > &)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
dictionary dict