35 template<
class CloudType>
44 Random&
rndGen = cloud.rndGen();
46 scalar ChiAMinusOne = ChiA - 1;
47 scalar ChiBMinusOne = ChiB - 1;
49 if (ChiAMinusOne < SMALL && ChiBMinusOne < SMALL)
51 return rndGen.sample01<scalar>();
61 energyRatio =
rndGen.sample01<scalar>();
63 if (ChiAMinusOne < SMALL)
65 P = 1.0 -
pow(energyRatio, ChiB);
67 else if (ChiBMinusOne < SMALL)
69 P = 1.0 -
pow(energyRatio, ChiA);
76 (ChiAMinusOne + ChiBMinusOne)*energyRatio/ChiAMinusOne,
81 (ChiAMinusOne + ChiBMinusOne)*(1 - energyRatio)
86 }
while (P <
rndGen.sample01<scalar>());
94 template<
class CloudType>
103 Tref_(this->coeffDict().getScalar(
"Tref")),
104 relaxationCollisionNumber_
106 this->coeffDict().getScalar(
"relaxationCollisionNumber")
113 template<
class CloudType>
121 template<
class CloudType>
128 template<
class CloudType>
137 label typeIdP = pP.typeId();
138 label typeIdQ = pQ.typeId();
143 cloud.constProps(typeIdP).d()
144 +
cloud.constProps(typeIdQ).d()
150 cloud.constProps(typeIdP).omega()
151 +
cloud.constProps(typeIdQ).omega()
154 scalar cR =
mag(pP.U() - pQ.U());
161 scalar mP = cloud.constProps(typeIdP).mass();
162 scalar mQ = cloud.constProps(typeIdQ).mass();
163 scalar mR = mP*mQ/(mP + mQ);
175 template<
class CloudType>
184 label typeIdP = pP.typeId();
185 label typeIdQ = pQ.typeId();
188 scalar& EiP = pP.Ei();
189 scalar& EiQ = pQ.Ei();
193 scalar inverseCollisionNumber = 1/relaxationCollisionNumber_;
199 scalar preCollisionEiP = EiP;
200 scalar preCollisionEiQ = EiQ;
202 direction iDofP =
cloud.constProps(typeIdP).internalDegreesOfFreedom();
203 direction iDofQ =
cloud.constProps(typeIdQ).internalDegreesOfFreedom();
208 cloud.constProps(typeIdP).omega()
209 +
cloud.constProps(typeIdQ).omega()
212 scalar mP =
cloud.constProps(typeIdP).mass();
213 scalar mQ =
cloud.constProps(typeIdQ).mass();
214 scalar mR = mP*mQ/(mP + mQ);
215 vector Ucm = (mP*UP + mQ*UQ)/(mP + mQ);
216 scalar cRsqr =
magSqr(UP - UQ);
217 scalar availableEnergy = 0.5*mR*cRsqr;
218 scalar ChiB = 2.5 - omegaPQ;
222 if (inverseCollisionNumber >
rndGen.sample01<scalar>())
224 availableEnergy += preCollisionEiP;
229 1.0 -
pow(
rndGen.sample01<scalar>(), (1.0/ChiB));
230 EiP = energyRatio*availableEnergy;
234 scalar ChiA = 0.5*iDofP;
235 EiP = energyRatio(ChiA, ChiB)*availableEnergy;
238 availableEnergy -= EiP;
244 if (inverseCollisionNumber >
rndGen.sample01<scalar>())
246 availableEnergy += preCollisionEiQ;
250 const scalar energyRatio =
251 1.0 -
pow(
rndGen.sample01<scalar>(), (1.0/ChiB));
253 EiQ = energyRatio*availableEnergy;
257 const scalar ChiA = 0.5*iDofQ;
258 EiQ = energyRatio(ChiA, ChiB)*availableEnergy;
261 availableEnergy -= EiQ;
266 scalar cR =
sqrt(2.0*availableEnergy/mR);
269 scalar cosTheta = 2.0*
rndGen.sample01<scalar>() - 1.0;
270 scalar sinTheta =
sqrt(1.0 - cosTheta*cosTheta);
273 vector postCollisionRelU =
282 UP = Ucm + postCollisionRelU*mQ/(mP + mQ);
283 UQ = Ucm - postCollisionRelU*mP/(mP + mQ);
virtual void collide(typename CloudType::parcelType &pP, typename CloudType::parcelType &pQ)
Apply collision.
DSMCCloud< dsmcParcel > CloudType
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
virtual scalar sigmaTcR(const typename CloudType::parcelType &pP, const typename CloudType::parcelType &pQ) const
Return the collision cross section * relative velocity product.
A list of keyword definitions, which are a keyword followed by a number of values (eg...
virtual bool active() const
Flag to indicate whether model activates collision model.
Templated DSMC particle collision class.
dimensionedScalar sqrt(const dimensionedScalar &ds)
virtual ~LarsenBorgnakkeVariableHardSphere()
Destructor.
dimensionedScalar cos(const dimensionedScalar &ds)
dimensionedScalar exp(const dimensionedScalar &ds)
constexpr scalar twoPi(2 *M_PI)
A cloud is a registry collection of lagrangian particles.
constexpr scalar pi(M_PI)
ParcelType parcelType
Type of parcel the cloud was instantiated for.
dimensionedScalar sin(const dimensionedScalar &ds)
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
dimensionedScalar lgamma(const dimensionedScalar &ds)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
const dimensionedScalar k
Boltzmann constant.
LarsenBorgnakkeVariableHardSphere(const dictionary &dict, CloudType &cloud)
Construct from dictionary.
Templated base class for dsmc cloud.
dimensioned< typename typeOfMag< Type >::type > magSqr(const dimensioned< Type > &dt)
Variable Hard Sphere BinaryCollision Model with Larsen Borgnakke internal energy redistribution. Based on the INELRS subroutine in Bird's DSMC0R.FOR.