38 template<
class ParticleType>
41 const polyBoundaryMesh& pbm = polyMesh_.boundaryMesh();
45 if (isA<cyclicAMIPolyPatch>(pbm[patchi]))
47 const cyclicAMIPolyPatch& cami =
48 refCast<const cyclicAMIPolyPatch>(pbm[
patchi]);
52 ok = ok && (cami.AMI().singlePatchProc() != -1);
60 <<
"Particle tracking across AMI patches is only currently " 61 <<
"supported for cases where the AMI patches reside on a " 67 template<
class ParticleType>
70 cellWallFacesPtr_.reset(
new PackedBoolList(pMesh().nCells(),
false));
72 PackedBoolList& cellWallFaces = cellWallFacesPtr_();
74 const polyBoundaryMesh& patches = polyMesh_.boundaryMesh();
78 if (isA<wallPolyPatch>(patches[patchi]))
80 const polyPatch& patch = patches[
patchi];
82 const labelList& pFaceCells = patch.faceCells();
86 cellWallFaces[pFaceCells[pFCI]] =
true;
95 template<
class ParticleType>
99 const word& cloudName,
103 cloud(pMesh, cloudName),
109 globalPositionsPtr_()
116 polyMesh_.tetBasePtIs();
118 if (particles.size())
127 template<
class ParticleType>
131 if (!cellWallFacesPtr_.valid())
136 return cellWallFacesPtr_();
140 template<
class ParticleType>
147 template<
class ParticleType>
150 delete(this->
remove(&
p));
154 template<
class ParticleType>
159 ParticleType&
p = pIter();
164 <<
"deleting lost particle at position " << p.position()
173 template<
class ParticleType>
178 ParticleType::particleCount_ = 0;
183 template<
class ParticleType>
184 template<
class TrackData>
201 const labelList& neighbourProcs = pData[Pstream::myProcNo()];
204 labelList neighbourProcIndices(Pstream::nProcs(), -1);
208 neighbourProcIndices[neighbourProcs[i]] = i;
214 pIter().stepFraction() = 0;
218 nTrackingRescues_ = 0;
225 neighbourProcs.
size()
232 neighbourProcs.
size()
239 globalPositionsPtr_.clear();
245 forAll(patchIndexTransferLists, i)
247 patchIndexTransferLists[i].clear();
253 ParticleType&
p = pIter();
256 bool keepParticle = p.move(td, trackTime);
267 && td.switchProcessor
275 if (procPatchIndices[patchi] != -1)
277 label n = neighbourProcIndices
279 refCast<const processorPolyPatch>
285 p.prepareForParallelTransfer(patchi, td);
287 particleTransferLists[
n].
append(this->
remove(&p));
289 patchIndexTransferLists[
n].append
291 procPatchNeighbours[patchi]
302 if (!Pstream::parRun())
312 forAll(particleTransferLists, i)
314 if (particleTransferLists[i].size())
323 << patchIndexTransferLists[i]
324 << particleTransferLists[i];
331 pBufs.finishedSends(allNTrans);
334 bool transfered =
false;
354 label neighbProci = neighbourProcs[i];
356 label nRec = allNTrans[neighbProci];
360 UIPstream particleStream(neighbProci, pBufs);
362 labelList receivePatchIndex(particleStream);
367 typename ParticleType::iNew(polyMesh_)
374 ParticleType& newp = newpIter();
376 label patchi = procPatches[receivePatchIndex[pI++]];
378 newp.correctAfterParallelTransfer(patchi, td);
380 addParticle(newParticles.
remove(&newp));
390 if (nTrackingRescues_ > 0)
392 Info<< nTrackingRescues_ <<
" tracking rescue corrections" <<
endl;
398 template<
class ParticleType>
401 if (!globalPositionsPtr_.valid())
404 <<
"Global positions are not available. " 405 <<
"Cloud::storeGlobalPositions has not been called." 411 cellWallFacesPtr_.clear();
416 polyMesh_.tetBasePtIs();
418 const vectorField& positions = globalPositionsPtr_();
423 iter().autoMap(positions[i], mapper);
429 template<
class ParticleType>
434 this->db().time().
path()/this->
name() +
"_positions.obj" 439 const ParticleType&
p = pIter();
440 pObj<<
"v " << p.position().x() <<
" " << p.position().y() <<
" " 441 << p.position().z() <<
nl;
448 template<
class ParticleType>
456 globalPositionsPtr_.reset(
new vectorField(this->size()));
463 positions[i] = iter().position();
Template class for intrusive linked lists.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
const labelList & processorPatchIndices() const
Return list of indices into processorPatches_ for each patch.
void cloudReset(const Cloud< ParticleType > &c)
Reset the particles.
#define forAll(list, i)
Loop across all elements in list.
const labelList & processorPatches() const
Return list of processor patch labels.
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
errorManipArg< error, int > exit(error &err, const int errNo=1)
void deleteParticle(ParticleType &)
Remove particle from cloud and delete.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
label nInternalFaces() const
Various mesh related information for a parallel run. Upon construction, constructs all info using par...
#define forAllIter(Container, container, iter)
Iterate across all elements in the container object of type.
void size(const label)
Override size to be inconsistent with allocated storage.
Cloud(const polyMesh &mesh, const word &cloudName, const IDLList< ParticleType > &particles)
Construct from mesh and a list of particles.
Ostream & endl(Ostream &os)
Add newline and flush stream.
volVectorField vectorField(fieldObject, mesh)
void autoMap(const mapPolyMesh &)
Remap the cells of particles corresponding to the.
Combination-Reduction operation for a parallel run. The information from all nodes is collected on th...
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Input inter-processor communications stream operating on external buffer.
void addParticle(ParticleType *pPtr)
Transfer particle to cloud.
void storeGlobalPositions() const
Call this before a topology change. Stores the particles global.
void clear()
Clear the list, i.e. set size to zero.
A class for handling words, derived from string.
A cloud is a collection of lagrangian particles.
void append(const T &)
Append an element at the end of the list.
virtual void flush()
Flush stream.
List< label > labelList
A List of labels.
T * remove(T *p)
Remove and return element.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
errorManip< error > abort(error &err)
Base cloud calls templated on particle type.
void move(TrackData &td, const scalar trackTime)
Move the particles.
Output inter-processor communications stream operating on external buffer.
const PackedBoolList & cellHasWallFaces() const
Whether each cell has any wall faces (demand driven data)
fileName path(UMean.rootPath()/UMean.caseName()/"graphs"/UMean.instance())
const labelList & processorPatchNeighbours() const
Return processorPatchIndices of the neighbours.
Buffers for inter-processor communications streams (UOPstream, UIPstream).
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
word name(const complex &)
Return a string representation of a complex.
void deleteLostParticles()
Remove lost particles from cloud and delete.
#define WarningInFunction
Report a warning using Foam::Warning.
const dimensionedScalar c
Speed of light in a vacuum.
Mesh consisting of general polyhedral cells.
label whichPatch(const label faceIndex) const
Return patch index for a given face label.
void writePositions() const
Write positions to <cloudName>_positions.obj file.