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]);
50 ok = ok && cami.singlePatchProc() != -1;
57 <<
"Particle tracking across AMI patches is only currently " 58 <<
"supported for cases where the AMI patches reside on a " 66 template<
class ParticleType>
70 const word& cloudName,
74 cloud(pMesh, cloudName),
84 polyMesh_.tetBasePtIs();
85 polyMesh_.oldCellCentres();
96 template<
class ParticleType>
103 template<
class ParticleType>
106 delete(this->
remove(&
p));
110 template<
class ParticleType>
115 ParticleType&
p = pIter();
120 <<
"deleting lost particle at position " << p.position()
129 template<
class ParticleType>
134 ParticleType::particleCount_ = 0;
139 template<
class ParticleType>
140 template<
class TrackCloudType>
143 TrackCloudType& cloud,
144 typename ParticleType::trackingData& td,
145 const scalar trackTime
159 const labelList& neighbourProcs = pData[Pstream::myProcNo()];
162 labelList neighbourProcIndices(Pstream::nProcs(), -1);
166 neighbourProcIndices[neighbourProcs[i]] = i;
179 neighbourProcs.
size()
186 neighbourProcs.
size()
193 globalPositionsPtr_.clear();
199 forAll(patchIndexTransferLists, i)
201 patchIndexTransferLists[i].clear();
207 ParticleType&
p = pIter();
210 bool keepParticle = p.move(cloud, td, trackTime);
216 if (td.switchProcessor)
222 || !p.onBoundaryFace()
223 || procPatchNeighbours[p.patch()] < 0
227 <<
"Switch processor flag is true when no parallel " 228 <<
"transfer is possible. This is a bug." 233 const label patchi = p.patch();
235 const label n = neighbourProcIndices
237 refCast<const processorPolyPatch>
243 p.prepareForParallelTransfer();
245 particleTransferLists[
n].
append(this->
remove(&p));
247 patchIndexTransferLists[
n].append
249 procPatchNeighbours[patchi]
259 if (!Pstream::parRun())
269 forAll(particleTransferLists, i)
271 if (particleTransferLists[i].size())
280 << patchIndexTransferLists[i]
281 << particleTransferLists[i];
288 pBufs.finishedSends(allNTrans);
291 bool transferred =
false;
311 label neighbProci = neighbourProcs[i];
313 label nRec = allNTrans[neighbProci];
317 UIPstream particleStream(neighbProci, pBufs);
319 labelList receivePatchIndex(particleStream);
324 typename ParticleType::iNew(polyMesh_)
331 ParticleType& newp = newpIter();
333 label patchi = procPatches[receivePatchIndex[pI++]];
335 newp.correctAfterParallelTransfer(patchi, td);
337 addParticle(newParticles.
remove(&newp));
345 template<
class ParticleType>
348 if (!globalPositionsPtr_.valid())
351 <<
"Global positions are not available. " 352 <<
"Cloud::storeGlobalPositions has not been called." 359 polyMesh_.tetBasePtIs();
360 polyMesh_.oldCellCentres();
362 const vectorField& positions = globalPositionsPtr_();
367 iter().autoMap(positions[i], mapper);
373 template<
class ParticleType>
378 this->db().time().
path()/this->
name() +
"_positions.obj" 383 const ParticleType&
p = pIter();
384 pObj<<
"v " << p.position().x() <<
" " << p.position().y() <<
" " 385 << p.position().z() <<
nl;
392 template<
class ParticleType>
400 globalPositionsPtr_.reset(
new vectorField(this->size()));
407 positions[i] = iter().position();
Template class for intrusive linked lists.
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
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.
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.
void append(const T &)
Append an element at the end of the list.
virtual void flush()
Flush stream.
T * remove(T *p)
Remove and return element.
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
errorManip< error > abort(error &err)
rAUs append(new volScalarField(IOobject::groupName("rAU", phase1.name()), 1.0/(U1Eqn.A()+byDt(max(phase1.residualAlpha() - alpha1, scalar(0)) *rho1))))
Base cloud calls templated on particle type.
Output inter-processor communications stream operating on external buffer.
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.
void move(TrackCloudType &cloud, typename ParticleType::trackingData &td, const scalar trackTime)
Move the particles.
#define WarningInFunction
Report a warning using Foam::Warning.
const dimensionedScalar c
Speed of light in a vacuum.
Mesh consisting of general polyhedral cells.
void writePositions() const
Write positions to <cloudName>_positions.obj file.
fileName path(UMean.rootPath()/UMean.caseName()/functionObjects::writeFile::outputPrefix/"graphs"/UMean.instance())