Go to the documentation of this file.
40 template<
class ParticleType>
43 const polyBoundaryMesh& pbm = polyMesh_.boundaryMesh();
45 for (
const polyPatch& pp : pbm)
47 if (isA<cyclicAMIPolyPatch>(pp))
49 const cyclicAMIPolyPatch& cami =
50 refCast<const cyclicAMIPolyPatch>(pp);
54 ok = ok && (cami.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)
91 polyMesh_.tetBasePtIs();
102 template<
class ParticleType>
109 template<
class ParticleType>
112 delete(this->remove(&
p));
116 template<
class ParticleType>
119 for (ParticleType&
p : *
this)
124 <<
"deleting lost particle at position " <<
p.position()
133 template<
class ParticleType>
138 ParticleType::particleCount_ = 0;
143 template<
class ParticleType>
144 template<
class TrackCloudType>
147 TrackCloudType&
cloud,
148 typename ParticleType::trackingData& td,
149 const scalar trackTime
163 const labelList& neighbourProcs = pData[Pstream::myProcNo()];
166 labelList neighbourProcIndices(Pstream::nProcs(), -1);
170 neighbourProcIndices[neighbourProcs[i]] = i;
176 pIter().stepFraction() = 0;
183 neighbourProcs.size()
190 neighbourProcs.size()
197 globalPositionsPtr_.clear();
203 forAll(patchIndexTransferLists, i)
205 patchIndexTransferLists[i].clear();
209 for (ParticleType&
p : *
this)
212 bool keepParticle =
p.move(
cloud, td, trackTime);
218 if (td.switchProcessor)
224 || !
p.onBoundaryFace()
225 || procPatchNeighbours[
p.patch()] < 0
229 <<
"Switch processor flag is true when no parallel "
230 <<
"transfer is possible. This is a bug."
235 const label patchi =
p.patch();
237 const label
n = neighbourProcIndices
239 refCast<const processorPolyPatch>
245 p.prepareForParallelTransfer();
247 particleTransferLists[
n].append(this->remove(&
p));
249 patchIndexTransferLists[
n].append
251 procPatchNeighbours[patchi]
261 if (!Pstream::parRun())
271 forAll(particleTransferLists, i)
273 if (particleTransferLists[i].size())
282 << patchIndexTransferLists[i]
283 << particleTransferLists[i];
290 pBufs.finishedSends(allNTrans);
293 bool transferred =
false;
295 for (
const label
n : allNTrans)
311 for (
const label neighbProci : neighbourProcs)
313 label nRec = allNTrans[neighbProci];
317 UIPstream particleStream(neighbProci, pBufs);
319 labelList receivePatchIndex(particleStream);
324 typename ParticleType::iNew(polyMesh_)
329 for (ParticleType& newp : newParticles)
331 label patchi = procPatches[receivePatchIndex[pI++]];
333 newp.correctAfterParallelTransfer(patchi, td);
335 addParticle(newParticles.remove(&newp));
343 template<
class ParticleType>
346 if (!globalPositionsPtr_.valid())
349 <<
"Global positions are not available. "
350 <<
"Cloud::storeGlobalPositions has not been called."
356 cellWallFacesPtr_.clear();
361 polyMesh_.tetBasePtIs();
363 const vectorField& positions = globalPositionsPtr_();
368 iter().autoMap(positions[i], mapper);
374 template<
class ParticleType>
379 this->db().time().
path()/this->
name() +
"_positions.obj"
384 const ParticleType&
p = pIter();
385 const point position(
p.position());
386 pObj<<
"v " << position.x() <<
" " << position.y() <<
" "
387 << position.z() <<
nl;
394 template<
class ParticleType>
402 globalPositionsPtr_.reset(
new vectorField(this->size()));
409 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.
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.
rAUs append(new volScalarField(IOobject::groupName("rAU", phase1.name()), 1.0/(U1Eqn.A()+byDt(max(phase1.residualAlpha() - alpha1, scalar(0)) *rho1))))
Input inter-processor communications stream operating on external buffer.
#define WarningInFunction
Report a warning using Foam::Warning.
virtual void flush()
Flush stream.