Go to the documentation of this file.
41 template<
class ParticleType>
44 const polyBoundaryMesh& pbm = polyMesh_.boundaryMesh();
46 for (
const polyPatch& pp : pbm)
48 if (isA<cyclicAMIPolyPatch>(pp))
50 const cyclicAMIPolyPatch& cami =
51 refCast<const cyclicAMIPolyPatch>(pp);
55 ok = ok && (cami.AMI().singlePatchProc() != -1);
63 <<
"Particle tracking across AMI patches is only currently "
64 <<
"supported for cases where the AMI patches reside on a "
72 template<
class ParticleType>
84 globalPositionsPtr_(),
85 geometryType_(cloud::geometryType::COORDINATES)
88 polyMesh_.oldCellCentres();
93 polyMesh_.tetBasePtIs();
104 template<
class ParticleType>
111 template<
class ParticleType>
114 delete(this->remove(&
p));
118 template<
class ParticleType>
121 for (ParticleType&
p : *
this)
126 <<
"deleting lost particle at position " <<
p.position()
135 template<
class ParticleType>
140 ParticleType::particleCount_ = 0;
145 template<
class ParticleType>
146 template<
class TrackCloudType>
149 TrackCloudType&
cloud,
150 typename ParticleType::trackingData& td,
151 const scalar trackTime
165 const labelList& neighbourProcs = pData[Pstream::myProcNo()];
168 labelList neighbourProcIndices(Pstream::nProcs(), -1);
172 neighbourProcIndices[neighbourProcs[i]] = i;
185 neighbourProcs.size()
192 neighbourProcs.size()
199 globalPositionsPtr_.clear();
205 forAll(patchIndexTransferLists, i)
207 patchIndexTransferLists[i].clear();
211 for (ParticleType&
p : *
this)
214 bool keepParticle =
p.move(
cloud, td, trackTime);
220 if (td.switchProcessor)
226 || !
p.onBoundaryFace()
227 || procPatchNeighbours[
p.patch()] < 0
231 <<
"Switch processor flag is true when no parallel "
232 <<
"transfer is possible. This is a bug."
237 const label patchi =
p.patch();
239 const label
n = neighbourProcIndices
241 refCast<const processorPolyPatch>
247 p.prepareForParallelTransfer();
249 particleTransferLists[
n].append(this->remove(&
p));
251 patchIndexTransferLists[
n].append
253 procPatchNeighbours[patchi]
263 if (!Pstream::parRun())
273 forAll(particleTransferLists, i)
275 if (particleTransferLists[i].size())
284 << patchIndexTransferLists[i]
285 << particleTransferLists[i];
292 pBufs.finishedSends(allNTrans);
295 bool transferred =
false;
297 for (
const label
n : allNTrans)
313 for (
const label neighbProci : neighbourProcs)
315 label nRec = allNTrans[neighbProci];
319 UIPstream particleStream(neighbProci, pBufs);
321 labelList receivePatchIndex(particleStream);
326 typename ParticleType::iNew(polyMesh_)
331 for (ParticleType& newp : newParticles)
333 label patchi = procPatches[receivePatchIndex[pI++]];
335 newp.correctAfterParallelTransfer(patchi, td);
337 addParticle(newParticles.remove(&newp));
345 template<
class ParticleType>
348 if (!globalPositionsPtr_)
351 <<
"Global positions are not available. "
352 <<
"Cloud::storeGlobalPositions has not been called."
358 cellWallFacesPtr_.clear();
363 polyMesh_.tetBasePtIs();
364 polyMesh_.oldCellCentres();
366 const vectorField& positions = globalPositionsPtr_();
371 iter().autoMap(positions[i], mapper);
377 template<
class ParticleType>
382 this->db().time().
path()/this->
name() +
"_positions.obj"
387 const ParticleType&
p = pIter();
388 const point position(
p.position());
389 pObj<<
"v " << position.x() <<
" " << position.y() <<
" "
390 << position.z() <<
nl;
397 template<
class ParticleType>
405 globalPositionsPtr_.reset(
new vectorField(this->size()));
412 positions[i] = iter().position();
void addParticle(ParticleType *pPtr)
Transfer particle to cloud.
const word cloudName(propsDict.get< word >("cloud"))
A class for handling words, derived from Foam::string.
Output inter-processor communications stream operating on external buffer.
A polyBoundaryMesh is a polyPatch list with additional search methods and registered IO.
const labelList & processorPatches() const
Return list of processor patch labels.
void writePositions() const
Write positions to <cloudName>_positions.obj file.
Buffers for inter-processor communications streams (UOPstream, UIPstream).
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Ostream & endl(Ostream &os)
Add newline and flush stream.
rAUs append(new volScalarField(IOobject::groupName("rAU", phase1.name()), 1.0/(U1Eqn.A()+byDt(max(phase1.residualAlpha() - alpha1, scalar(0)) *rho1))))
Mesh consisting of general polyhedral cells.
#define forAll(list, i)
Loop across all elements in list.
Field< vector > vectorField
Specialisation of Field<T> for vector.
Combination-Reduction operation for a parallel run. The information from all nodes is collected on th...
word name(const complex &c)
Return string representation of complex.
void storeGlobalPositions() const
Call this before a topology change.
void deleteParticle(ParticleType &p)
Remove particle from cloud and delete.
void autoMap(const mapPolyMesh &)
Remap the cells of particles corresponding to the.
void deleteLostParticles()
Remove lost particles from cloud and delete.
#define forAllIters(container, iter)
Iterate across all elements in the container object.
Template class for intrusive linked lists.
void move(TrackCloudType &cloud, typename ParticleType::trackingData &td, const scalar trackTime)
Move the particles.
Various mesh related information for a parallel run. Upon construction, constructs all info using par...
reduce(hasMovingMesh, orOp< bool >())
errorManip< error > abort(error &err)
const labelList & processorPatchNeighbours() const
Return processorPatchIndices of the neighbours.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Output to file stream, using an OSstream.
void cloudReset(const Cloud< ParticleType > &c)
Reset the particles.
A cloud is a registry collection of lagrangian particles.
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
forAllConstIters(mixture.phases(), phase)
Base cloud calls templated on particle type.
void clear()
Clear the list, i.e. set size to zero.
fileName path(UMean.rootPath()/UMean.caseName()/"graphs"/UMean.instance())
Cloud(const polyMesh &mesh, const word &cloudName, const IDLList< ParticleType > &particles)
Construct from mesh and a list of particles.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
const dimensionedScalar c
Speed of light in a vacuum.
Input inter-processor communications stream operating on external buffer.
#define WarningInFunction
Report a warning using Foam::Warning.
virtual void flush()
Flush stream.