32 template<
class CloudType>
35 template<
class CloudType>
42 template<
class CloudType>
48 typename CloudType::parcelType& p = iter();
57 template<
class CloudType>
60 PstreamBuffers pBufs(Pstream::commsTypes::nonBlocking);
62 label startOfRequests = Pstream::nRequests();
66 realRealInteraction();
68 il_.receiveReferredData(pBufs, startOfRequests);
70 realReferredInteraction();
74 template<
class CloudType>
80 typename CloudType::parcelType* pA_ptr =
nullptr;
81 typename CloudType::parcelType* pB_ptr =
nullptr;
83 List<DynamicList<typename CloudType::parcelType*>>& cellOccupancy =
84 this->owner().cellOccupancy();
89 forAll(cellOccupancy[realCelli], a)
91 pA_ptr = cellOccupancy[realCelli][a];
93 forAll(dil[realCelli], interactingCells)
95 List<typename CloudType::parcelType*> cellBParcels =
96 cellOccupancy[dil[realCelli][interactingCells]];
101 pB_ptr = cellBParcels[
b];
103 evaluatePair(*pA_ptr, *pB_ptr);
108 forAll(cellOccupancy[realCelli], aO)
110 pB_ptr = cellOccupancy[realCelli][aO];
116 evaluatePair(*pA_ptr, *pB_ptr);
124 template<
class CloudType>
130 List<IDLList<typename CloudType::parcelType>>& referredParticles =
131 il_.referredParticles();
133 List<DynamicList<typename CloudType::parcelType*>>& cellOccupancy =
134 this->owner().cellOccupancy();
139 IDLList<typename CloudType::parcelType>& refCellRefParticles =
140 referredParticles[refCelli];
142 const labelList& realCells = ril[refCelli];
148 typename IDLList<typename CloudType::parcelType>,
156 forAll(realCells, realCelli)
158 List<typename CloudType::parcelType*> realCellParcels =
159 cellOccupancy[realCells[realCelli]];
161 forAll(realCellParcels, realParcelI)
165 *realCellParcels[realParcelI],
175 template<
class CloudType>
178 const polyMesh& mesh = this->owner().mesh();
184 const labelList& patchID = mesh.boundaryMesh().patchID();
188 List<DynamicList<typename CloudType::parcelType*>>& cellOccupancy =
189 this->owner().cellOccupancy();
192 DynamicList<point> flatSitePoints;
193 DynamicList<scalar> flatSiteExclusionDistancesSqr;
194 DynamicList<WallSiteData<vector>> flatSiteData;
195 DynamicList<point> otherSitePoints;
196 DynamicList<scalar> otherSiteDistances;
197 DynamicList<WallSiteData<vector>> otherSiteData;
198 DynamicList<point> sharpSitePoints;
199 DynamicList<scalar> sharpSiteExclusionDistancesSqr;
200 DynamicList<WallSiteData<vector>> sharpSiteData;
205 const labelList& realWallFaces = directWallFaces[realCelli];
208 forAll(cellOccupancy[realCelli], cellParticleI)
210 flatSitePoints.clear();
211 flatSiteExclusionDistancesSqr.clear();
212 flatSiteData.clear();
213 otherSitePoints.clear();
214 otherSiteDistances.clear();
215 otherSiteData.clear();
216 sharpSitePoints.clear();
217 sharpSiteExclusionDistancesSqr.clear();
218 sharpSiteData.clear();
220 typename CloudType::parcelType& p =
221 *cellOccupancy[realCelli][cellParticleI];
223 const point& pos = p.position();
225 scalar r = wallModel_->pREff(p);
229 forAll(realWallFaces, realWallFacei)
231 label realFacei = realWallFaces[realWallFacei];
233 pointHit nearest = mesh.faces()[realFacei].nearestPoint
239 if (nearest.distance() < r)
241 vector normal = mesh.faceAreas()[realFacei];
243 normal /=
mag(normal);
245 const vector& nearPt = nearest.rawPoint();
249 scalar normalAlignment = normal & pW/(
mag(pW) + small);
252 label patchi = patchID[realFacei - mesh.nInternalFaces()];
255 realFacei - mesh.boundaryMesh()[
patchi].start();
257 WallSiteData<vector> wSD
260 U.boundaryField()[
patchi][patchFacei]
263 if (normalAlignment > cosPhiMinFlatWall)
271 !duplicatePointInList
275 sqr(r*flatWallDuplicateExclusion)
279 flatSitePoints.append(nearPt);
281 flatSiteExclusionDistancesSqr.append
283 sqr(r) -
sqr(nearest.distance())
286 flatSiteData.append(wSD);
291 otherSitePoints.append(nearPt);
293 otherSiteDistances.append(nearest.distance());
295 otherSiteData.append(wSD);
303 const labelList& cellRefWallFaces = il_.rwfilInverse()[realCelli];
305 forAll(cellRefWallFaces, rWFI)
307 label refWallFacei = cellRefWallFaces[rWFI];
309 const referredWallFace& rwf =
310 il_.referredWallFaces()[refWallFacei];
314 pointHit nearest = rwf.nearestPoint(pos, pts);
316 if (nearest.distance() < r)
318 const vector normal = rwf.normal(pts);
319 const vector& nearPt = nearest.rawPoint();
323 scalar normalAlignment = normal & pW/
mag(pW);
327 WallSiteData<vector> wSD
330 il_.referredWallData()[refWallFacei]
333 if (normalAlignment > cosPhiMinFlatWall)
341 !duplicatePointInList
345 sqr(r*flatWallDuplicateExclusion)
349 flatSitePoints.append(nearPt);
351 flatSiteExclusionDistancesSqr.append
353 sqr(r) -
sqr(nearest.distance())
356 flatSiteData.append(wSD);
361 otherSitePoints.append(nearPt);
363 otherSiteDistances.append(nearest.distance());
365 otherSiteData.append(wSD);
381 sortedOrder(otherSiteDistances, sortedOtherSiteIndices);
383 forAll(sortedOtherSiteIndices, siteI)
385 label orderedIndex = sortedOtherSiteIndices[siteI];
387 const point& otherPt = otherSitePoints[orderedIndex];
391 !duplicatePointInList
395 flatSiteExclusionDistancesSqr
404 !duplicatePointInList
408 sharpSiteExclusionDistancesSqr
412 sharpSitePoints.append(otherPt);
414 sharpSiteExclusionDistancesSqr.append
416 sqr(r) -
sqr(otherSiteDistances[orderedIndex])
419 sharpSiteData.append(otherSiteData[orderedIndex]);
437 template<
class CloudType>
440 const DynamicList<point>& existingPoints,
441 const point& pointToTest,
442 scalar duplicateRangeSqr
447 if (
magSqr(existingPoints[i] - pointToTest) < duplicateRangeSqr)
457 template<
class CloudType>
460 const DynamicList<point>& existingPoints,
461 const point& pointToTest,
467 if (
magSqr(existingPoints[i] - pointToTest) < duplicateRangeSqr[i])
477 template<
class CloudType>
484 typename CloudType::parcelType& p = iter();
486 p.collisionRecords().update();
491 template<
class CloudType>
494 typename CloudType::parcelType& pA,
495 typename CloudType::parcelType& pB
498 pairModel_->evaluatePair(pA, pB);
502 template<
class CloudType>
505 typename CloudType::parcelType& p,
506 const List<point>& flatSitePoints,
507 const List<WallSiteData<vector>>& flatSiteData,
508 const List<point>& sharpSitePoints,
509 const List<WallSiteData<vector>>& sharpSiteData
512 wallModel_->evaluateWall
525 template<
class CloudType>
552 this->coeffDict().template lookup<scalar>(
"maxInteractionDistance"),
555 this->coeffDict().lookupOrDefault
557 "writeReferredParticleCloud",
561 this->coeffDict().lookupOrDefault(
"U",
word(
"U"))
566 template<
class CloudType>
575 il_(cm.
owner().mesh())
584 template<
class CloudType>
591 template<
class CloudType>
594 label nSubCycles = 1;
596 if (pairModel_->controlsTimestep())
603 nSubCycles =
max(nSubCycles, nPairSubCycles);
606 if (wallModel_->controlsTimestep())
613 nSubCycles =
max(nSubCycles, nWallSubCycles);
620 template<
class CloudType>
List< labelList > labelListList
A List of labelList.
#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.
void sortedOrder(const UList< T > &, labelList &order)
Generate the (stable) sort order for the list.
DSMCCloud< dsmcParcel > CloudType
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)
A simple wrapper around bool so that it can be read as a word: true/false, on/off, yes/no, y/n, t/f, or none/any.
Vector< scalar > vector
A scalar version of the templated Vector.
GeometricField< vector, fvPatchField, volMesh > volVectorField
const CloudType & owner() const
Return const access to the owner cloud.
dimensionedScalar pos(const dimensionedScalar &ds)
vectorField pointField
pointField is a vectorField.
PairCollision(const dictionary &dict, CloudType &owner)
Construct from components.
virtual label nSubCycles() const
Return the number of times to subcycle the current.
const dimensionedScalar & b
Wien displacement law constant: default SI units: [m K].
A class for handling words, derived from string.
List< scalar > scalarList
A List of scalars.
Templated pair interaction class.
Templated collision model class.
List< label > labelList
A List of labels.
Templated wall interaction class.
dimensioned< scalar > magSqr(const dimensioned< Type > &)
const List< DynamicList< molecule * > > & cellOccupancy
vector point
Point is a vector.
dimensioned< scalar > mag(const dimensioned< Type > &)
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
PointHit< point > pointHit
virtual ~PairCollision()
Destructor.
#define NotImplemented
Issue a FatalErrorIn for a function not currently implemented.
Templated base class for dsmc cloud.