35 template<
class CloudType>
38 typename CloudType::parcelType::trackingData& td,
46 occupancy[iter().cell()]++;
58 pInCell(iter().
cell(), occupancy[iter().
cell()]++) = &iter();
61 for (
label celli=0; celli<this->owner().mesh().nCells(); celli++)
65 if (pInCelli.
size() >= 2)
74 scalar m1 = p1.nParticle()*p1.mass();
75 scalar m2 = p2.nParticle()*p2.mass();
77 bool massChanged = collideParcels(dt, p1, p2, m1, m2);
84 p1.setCellValues(this->owner(), td);
85 p1.rho() = liquids_.rho(td.pc(), p1.T(), X);
86 p1.Cp() = liquids_.Cp(td.pc(), p1.T(), X);
87 p1.sigma() = liquids_.sigma(td.pc(), p1.T(), X);
88 p1.mu() = liquids_.mu(td.pc(), p1.T(), X);
89 p1.d() =
cbrt(6.0*m1/(p1.nParticle()*p1.rho()*
pi));
95 p2.setCellValues(this->owner(), td);
96 p2.rho() = liquids_.rho(td.pc(), p2.T(), X);
97 p2.Cp() = liquids_.Cp(td.pc(), p2.T(), X);
98 p2.sigma() = liquids_.sigma(td.pc(), p2.T(), X);
99 p2.mu() = liquids_.mu(td.pc(), p2.T(), X);
100 p2.d() =
cbrt(6.0*m2/(p2.nParticle()*p2.rho()*
pi));
112 scalar mass = p.nParticle()*p.mass();
114 if (mass < this->owner().constProps().minParcelMass())
116 this->owner().deleteParticle(p);
122 template<
class CloudType>
133 if ((m1 < rootVSmall) || (m2 < rootVSmall))
138 const scalar Vc = this->owner().mesh().V()[p1.cell()];
139 const scalar d1 = p1.d();
140 const scalar d2 = p2.d();
142 scalar magUrel =
mag(p1.U() - p2.U());
143 scalar sumD = d1 + d2;
145 scalar nMin =
min(p1.nParticle(), p2.nParticle());
146 scalar
nu = nMin*nu0;
147 scalar collProb =
exp(-nu);
148 scalar xx = this->owner().rndGen().template sample01<scalar>();
155 return collideSorted(dt, p1, p2, m1, m2);
159 return collideSorted(dt, p2, p1, m2, m1);
169 template<
class CloudType>
179 const scalar nP1 = p1.nParticle();
180 const scalar nP2 = p2.nParticle();
182 const scalar sigma1 = p1.sigma();
183 const scalar sigma2 = p2.sigma();
185 const scalar d1 = p1.d();
186 const scalar d2 = p2.d();
188 const scalar T1 = p1.T();
189 const scalar T2 = p2.T();
191 const scalar
rho1 = p1.rho();
192 const scalar
rho2 = p2.rho();
199 scalar magURel =
mag(URel);
201 scalar
mTot = m1 + m2;
203 scalar gamma = d1/
max(rootVSmall, d2);
204 scalar
f =
pow3(gamma) + 2.7*gamma - 2.4*
sqr(gamma);
207 scalar Tave = (T1*m1 + T2*m2)/mTot;
210 scalar sigmaAve = sigma1;
211 if (
mag(T2 - T1) > small)
213 sigmaAve += (sigma2 - sigma1)*(Tave - T1)/(T2 - T1);
216 scalar Vtot = m1/rho1 + m2/
rho2;
217 scalar rhoAve = mTot/Vtot;
219 scalar dAve =
sqrt(d1*d2);
220 scalar WeColl = 0.5*rhoAve*
sqr(magURel)*dAve/
max(rootVSmall, sigmaAve);
222 scalar coalesceProb =
min(1.0, 2.4*f/
max(rootVSmall, WeColl));
224 scalar prob = this->owner().rndGen().template sample01<scalar>();
227 if (coalescence_ && prob < coalesceProb)
230 scalar nProb = prob*nP2/nP1;
236 scalar dm = nP1*nProb*m2/scalar(nP2);
241 p1.T() = (Tave*mTot - m2*T2)/m1;
243 p1.U() = (m1*U1 + (1.0 - m2/m2Org)*m2*U2)/m1;
245 p1.Y() = (m1Org*p1.Y() + dm*p2.Y())/m1;
247 p2.nParticle() = m2/(rho2*p2.volume());
254 scalar gf =
sqrt(prob) -
sqrt(coalesceProb);
255 scalar denom = 1.0 -
sqrt(coalesceProb);
270 vector v1p = (mr + m2*gf*URel)/mTot;
271 vector v2p = (mr - m1*gf*URel)/mTot;
276 p2.U() = (nP1*v2p + (nP2 - nP1)*U2)/nP2;
280 p1.U() = (nP2*v1p + (nP1 - nP2)*U1)/nP1;
291 template<
class CloudType>
296 const word& modelName
302 owner.db().template lookupObject<SLGThermo>(
"SLGThermo").liquids()
304 coalescence_(this->coeffDict().lookup(
"coalescence"))
308 template<
class CloudType>
322 template<
class CloudType>
ORourkeCollision(const dictionary &dict, CloudType &cloud, const word &modelName=typeName)
Construct from dictionary.
#define forAll(list, i)
Loop across all elements in list.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
A list of keyword definitions, which are a keyword followed by any number of values (e...
dimensioned< Type > max(const dimensioned< Type > &, const dimensioned< Type > &)
dimensionedSymmTensor sqr(const dimensionedVector &dv)
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
dimensionedScalar sqrt(const dimensionedScalar &ds)
virtual void collide(typename CloudType::parcelType::trackingData &td, const scalar dt)
Main collision routine.
virtual ~ORourkeCollision()
Destructor.
dimensionedScalar exp(const dimensionedScalar &ds)
A class for handling words, derived from string.
dimensionedScalar cbrt(const dimensionedScalar &ds)
virtual bool collideParcels(const scalar dt, parcelType &p1, parcelType &p2, scalar &m1, scalar &m2)
Collide parcels and return true if mass has changed.
A 1D vector of objects of type <T>, where the size of the vector is known and can be used for subscri...
A packed storage unstructured matrix of objects of type <T> using an offset table for access...
CloudType::parcelType parcelType
Convenience typedef to the cloud's parcel type.
virtual bool collideSorted(const scalar dt, parcelType &p1, parcelType &p2, scalar &m1, scalar &m2)
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
const liquidMixtureProperties & liquids_
Collision model by P.J. O'Rourke.
dimensionedScalar pow3(const dimensionedScalar &ds)
A cell is defined as a list of faces with extra functionality.
Templated stochastic collision model class.
dimensioned< scalar > mag(const dimensioned< Type > &)
const dimensionedScalar & rho2
label size() const
Return the number of elements in the UList.
const dimensionedScalar & rho1
Templated base class for dsmc cloud.
Switch coalescence_
Coalescence activation switch.