Go to the documentation of this file.
41 template<
class ParticleType>
45 for (
const polyPatch& pp : polyMesh_.boundaryMesh())
47 const auto* camipp = isA<cyclicAMIPolyPatch>(pp);
49 if (camipp && camipp->owner())
51 ok = (camipp->AMI().singlePatchProc() != -1);
62 <<
"Particle tracking across AMI patches is only currently "
63 <<
"supported for cases where the AMI patches reside on a "
71 template<
class ParticleType>
83 globalPositionsPtr_(),
84 geometryType_(cloud::geometryType::COORDINATES)
87 polyMesh_.oldCellCentres();
92 polyMesh_.tetBasePtIs();
103 template<
class ParticleType>
110 template<
class ParticleType>
113 delete(this->remove(&
p));
117 template<
class ParticleType>
120 for (ParticleType&
p : *
this)
125 <<
"deleting lost particle at position " <<
p.position()
134 template<
class ParticleType>
139 ParticleType::particleCount_ = 0;
144 template<
class ParticleType>
145 template<
class TrackCloudType>
148 TrackCloudType&
cloud,
149 typename ParticleType::trackingData& td,
150 const scalar trackTime
164 const labelList& neighbourProcs = pData[Pstream::myProcNo()];
167 labelList neighbourProcIndices(Pstream::nProcs(), -1);
171 neighbourProcIndices[neighbourProcs[i]] = i;
184 neighbourProcs.size()
191 neighbourProcs.size()
198 globalPositionsPtr_.clear();
204 forAll(patchIndexTransferLists, i)
206 patchIndexTransferLists[i].clear();
210 for (ParticleType&
p : *
this)
213 bool keepParticle =
p.move(
cloud, td, trackTime);
219 if (td.switchProcessor)
225 || !
p.onBoundaryFace()
226 || procPatchNeighbours[
p.patch()] < 0
230 <<
"Switch processor flag is true when no parallel "
231 <<
"transfer is possible. This is a bug."
236 const label patchi =
p.patch();
238 const label
n = neighbourProcIndices
240 refCast<const processorPolyPatch>
246 p.prepareForParallelTransfer();
248 particleTransferLists[
n].append(this->remove(&
p));
250 patchIndexTransferLists[
n].append
252 procPatchNeighbours[patchi]
262 if (!Pstream::parRun())
272 forAll(particleTransferLists, i)
274 if (particleTransferLists[i].size())
283 << patchIndexTransferLists[i]
284 << particleTransferLists[i];
291 pBufs.finishedSends(allNTrans);
294 bool transferred =
false;
296 for (
const label
n : allNTrans)
312 for (
const label neighbProci : neighbourProcs)
314 label nRec = allNTrans[neighbProci];
318 UIPstream particleStream(neighbProci, pBufs);
320 labelList receivePatchIndex(particleStream);
325 typename ParticleType::iNew(polyMesh_)
330 for (ParticleType& newp : newParticles)
332 label patchi = procPatches[receivePatchIndex[pI++]];
334 newp.correctAfterParallelTransfer(patchi, td);
336 addParticle(newParticles.remove(&newp));
344 template<
class ParticleType>
347 if (!globalPositionsPtr_)
350 <<
"Global positions are not available. "
351 <<
"Cloud::storeGlobalPositions has not been called."
357 cellWallFacesPtr_.clear();
362 polyMesh_.tetBasePtIs();
363 polyMesh_.oldCellCentres();
365 const vectorField& positions = globalPositionsPtr_();
370 iter().autoMap(positions[i], mapper);
376 template<
class ParticleType>
381 this->db().time().
path()/this->
name() +
"_positions.obj"
386 const ParticleType&
p = pIter();
387 const point position(
p.position());
388 pObj<<
"v " << position.x() <<
" " << position.y() <<
" "
389 << position.z() <<
nl;
396 template<
class ParticleType>
404 globalPositionsPtr_.reset(
new vectorField(this->size()));
411 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.
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.
const labelList & processorPatches() const noexcept
Return list of processor patch labels.
#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...
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)
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.
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Input inter-processor communications stream operating on external buffer.
#define WarningInFunction
Report a warning using Foam::Warning.
const labelList & processorPatchNeighbours() const noexcept
Return processorPatchIndices of the neighbours.
virtual void flush()
Flush stream.