Cloud.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2019 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "Cloud.H"
27 #include "processorPolyPatch.H"
28 #include "globalMeshData.H"
30 #include "mapPolyMesh.H"
31 #include "Time.H"
32 #include "OFstream.H"
33 #include "wallPolyPatch.H"
34 #include "cyclicAMIPolyPatch.H"
35 
36 // * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
37 
38 template<class ParticleType>
40 {
41  const polyBoundaryMesh& pbm = polyMesh_.boundaryMesh();
42  bool ok = true;
43  forAll(pbm, patchi)
44  {
45  if (isA<cyclicAMIPolyPatch>(pbm[patchi]))
46  {
47  const cyclicAMIPolyPatch& cami =
48  refCast<const cyclicAMIPolyPatch>(pbm[patchi]);
49 
50  ok = ok && cami.singlePatchProc() != -1;
51  }
52  }
53 
54  if (!ok)
55  {
57  << "Particle tracking across AMI patches is only currently "
58  << "supported for cases where the AMI patches reside on a "
59  << "single processor" << abort(FatalError);
60  }
61 }
62 
63 
64 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
65 
66 template<class ParticleType>
68 (
69  const polyMesh& pMesh,
70  const word& cloudName,
71  const IDLList<ParticleType>& particles
72 )
73 :
74  cloud(pMesh, cloudName),
76  polyMesh_(pMesh),
77  globalPositionsPtr_()
78 {
79  checkPatches();
80 
81  // Ask for the tetBasePtIs and oldCellCentres to trigger all processors to
82  // build them, otherwise, if some processors have no particles then there
83  // is a comms mismatch.
84  polyMesh_.tetBasePtIs();
85  polyMesh_.oldCellCentres();
86 
87  if (particles.size())
88  {
90  }
91 }
92 
93 
94 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
95 
96 template<class ParticleType>
98 {
99  this->append(pPtr);
100 }
101 
102 
103 template<class ParticleType>
105 {
106  delete(this->remove(&p));
107 }
108 
109 
110 template<class ParticleType>
112 {
113  forAllIter(typename Cloud<ParticleType>, *this, pIter)
114  {
115  ParticleType& p = pIter();
116 
117  if (p.cell() == -1)
118  {
120  << "deleting lost particle at position " << p.position()
121  << endl;
122 
123  deleteParticle(p);
124  }
125  }
126 }
127 
128 
129 template<class ParticleType>
131 {
132  // Reset particle count and particles only
133  // - not changing the cloud object registry or reference to the polyMesh
134  ParticleType::particleCount_ = 0;
136 }
137 
138 
139 template<class ParticleType>
140 template<class TrackCloudType>
142 (
143  TrackCloudType& cloud,
144  typename ParticleType::trackingData& td,
145  const scalar trackTime
146 )
147 {
148  const polyBoundaryMesh& pbm = pMesh().boundaryMesh();
149  const globalMeshData& pData = polyMesh_.globalData();
150 
151  // Which patches are processor patches
152  const labelList& procPatches = pData.processorPatches();
153 
154  // Indexing of equivalent patch on neighbour processor into the
155  // procPatches list on the neighbour
156  const labelList& procPatchNeighbours = pData.processorPatchNeighbours();
157 
158  // Which processors this processor is connected to
159  const labelList& neighbourProcs = pData[Pstream::myProcNo()];
160 
161  // Indexing from the processor number into the neighbourProcs list
162  labelList neighbourProcIndices(Pstream::nProcs(), -1);
163 
164  forAll(neighbourProcs, i)
165  {
166  neighbourProcIndices[neighbourProcs[i]] = i;
167  }
168 
169  // Initialise the stepFraction moved for the particles
170  forAllIter(typename Cloud<ParticleType>, *this, pIter)
171  {
172  pIter().reset();
173  }
174 
175  // List of lists of particles to be transferred for all of the
176  // neighbour processors
177  List<IDLList<ParticleType>> particleTransferLists
178  (
179  neighbourProcs.size()
180  );
181 
182  // List of destination processorPatches indices for all of the
183  // neighbour processors
184  List<DynamicList<label>> patchIndexTransferLists
185  (
186  neighbourProcs.size()
187  );
188 
189  // Allocate transfer buffers
190  PstreamBuffers pBufs(Pstream::commsTypes::nonBlocking);
191 
192  // Clear the global positions as there are about to change
193  globalPositionsPtr_.clear();
194 
195  // While there are particles to transfer
196  while (true)
197  {
198  particleTransferLists = IDLList<ParticleType>();
199  forAll(patchIndexTransferLists, i)
200  {
201  patchIndexTransferLists[i].clear();
202  }
203 
204  // Loop over all particles
205  forAllIter(typename Cloud<ParticleType>, *this, pIter)
206  {
207  ParticleType& p = pIter();
208 
209  // Move the particle
210  bool keepParticle = p.move(cloud, td, trackTime);
211 
212  // If the particle is to be kept
213  // (i.e. it hasn't passed through an inlet or outlet)
214  if (keepParticle)
215  {
216  if (td.switchProcessor)
217  {
218  #ifdef FULLDEBUG
219  if
220  (
221  !Pstream::parRun()
222  || !p.onBoundaryFace()
223  || procPatchNeighbours[p.patch()] < 0
224  )
225  {
227  << "Switch processor flag is true when no parallel "
228  << "transfer is possible. This is a bug."
229  << exit(FatalError);
230  }
231  #endif
232 
233  const label patchi = p.patch();
234 
235  const label n = neighbourProcIndices
236  [
237  refCast<const processorPolyPatch>
238  (
239  pbm[patchi]
240  ).neighbProcNo()
241  ];
242 
243  p.prepareForParallelTransfer();
244 
245  particleTransferLists[n].append(this->remove(&p));
246 
247  patchIndexTransferLists[n].append
248  (
249  procPatchNeighbours[patchi]
250  );
251  }
252  }
253  else
254  {
255  deleteParticle(p);
256  }
257  }
258 
259  if (!Pstream::parRun())
260  {
261  break;
262  }
263 
264 
265  // Clear transfer buffers
266  pBufs.clear();
267 
268  // Stream into send buffers
269  forAll(particleTransferLists, i)
270  {
271  if (particleTransferLists[i].size())
272  {
273  UOPstream particleStream
274  (
275  neighbourProcs[i],
276  pBufs
277  );
278 
279  particleStream
280  << patchIndexTransferLists[i]
281  << particleTransferLists[i];
282  }
283  }
284 
285 
286  // Start sending. Sets number of bytes transferred
287  labelList allNTrans(Pstream::nProcs());
288  pBufs.finishedSends(allNTrans);
289 
290 
291  bool transferred = false;
292 
293  forAll(allNTrans, i)
294  {
295  if (allNTrans[i])
296  {
297  transferred = true;
298  break;
299  }
300  }
301  reduce(transferred, orOp<bool>());
302 
303  if (!transferred)
304  {
305  break;
306  }
307 
308  // Retrieve from receive buffers
309  forAll(neighbourProcs, i)
310  {
311  label neighbProci = neighbourProcs[i];
312 
313  label nRec = allNTrans[neighbProci];
314 
315  if (nRec)
316  {
317  UIPstream particleStream(neighbProci, pBufs);
318 
319  labelList receivePatchIndex(particleStream);
320 
321  IDLList<ParticleType> newParticles
322  (
323  particleStream,
324  typename ParticleType::iNew(polyMesh_)
325  );
326 
327  label pI = 0;
328 
329  forAllIter(typename Cloud<ParticleType>, newParticles, newpIter)
330  {
331  ParticleType& newp = newpIter();
332 
333  label patchi = procPatches[receivePatchIndex[pI++]];
334 
335  newp.correctAfterParallelTransfer(patchi, td);
336 
337  addParticle(newParticles.remove(&newp));
338  }
339  }
340  }
341  }
342 }
343 
344 
345 template<class ParticleType>
347 {
348  if (!globalPositionsPtr_.valid())
349  {
351  << "Global positions are not available. "
352  << "Cloud::storeGlobalPositions has not been called."
353  << exit(FatalError);
354  }
355 
356  // Ask for the tetBasePtIs to trigger all processors to build
357  // them, otherwise, if some processors have no particles then
358  // there is a comms mismatch.
359  polyMesh_.tetBasePtIs();
360  polyMesh_.oldCellCentres();
361 
362  const vectorField& positions = globalPositionsPtr_();
363 
364  label i = 0;
365  forAllIter(typename Cloud<ParticleType>, *this, iter)
366  {
367  iter().autoMap(positions[i], mapper);
368  ++ i;
369  }
370 }
371 
372 
373 template<class ParticleType>
375 {
376  OFstream pObj
377  (
378  this->db().time().path()/this->name() + "_positions.obj"
379  );
380 
381  forAllConstIter(typename Cloud<ParticleType>, *this, pIter)
382  {
383  const ParticleType& p = pIter();
384  pObj<< "v " << p.position().x() << " " << p.position().y() << " "
385  << p.position().z() << nl;
386  }
387 
388  pObj.flush();
389 }
390 
391 
392 template<class ParticleType>
394 {
395  // Store the global positions for later use by autoMap. It would be
396  // preferable not to need this. If the mapPolyMesh object passed to autoMap
397  // had a copy of the old mesh then the global positions could be recovered
398  // within autoMap, and this pre-processing would not be necessary.
399 
400  globalPositionsPtr_.reset(new vectorField(this->size()));
401 
402  vectorField& positions = globalPositionsPtr_();
403 
404  label i = 0;
405  forAllConstIter(typename Cloud<ParticleType>, *this, iter)
406  {
407  positions[i] = iter().position();
408  ++ i;
409  }
410 }
411 
412 
413 // * * * * * * * * * * * * * * * * IOStream operators * * * * * * * * * * * //
414 
415 #include "CloudIO.C"
416 
417 // ************************************************************************* //
Template class for intrusive linked lists.
Definition: ILList.H:50
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:434
void cloudReset(const Cloud< ParticleType > &c)
Reset the particles.
Definition: Cloud.C:130
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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.
Definition: label.H:59
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
void deleteParticle(ParticleType &)
Remove particle from cloud and delete.
Definition: Cloud.C:104
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
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.
Definition: UList.H:459
Output to file stream.
Definition: OFstream.H:82
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
Cloud(const polyMesh &mesh, const word &cloudName, const IDLList< ParticleType > &particles)
Construct from mesh and a list of particles.
Definition: Cloud.C:68
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:256
volVectorField vectorField(fieldObject, mesh)
void autoMap(const mapPolyMesh &)
Remap the cells of particles corresponding to the.
Definition: Cloud.C:346
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.
Definition: mapPolyMesh.H:158
Input inter-processor communications stream operating on external buffer.
Definition: UIPstream.H:53
void addParticle(ParticleType *pPtr)
Transfer particle to cloud.
Definition: Cloud.C:97
void storeGlobalPositions() const
Call this before a topology change. Stores the particles global.
Definition: Cloud.C:393
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:124
A class for handling words, derived from string.
Definition: word.H:59
void append(const T &)
Append an element at the end of the list.
Definition: ListI.H:177
virtual void flush()
Flush stream.
Definition: OSstream.C:254
T * remove(T *p)
Remove and return element.
Definition: UILList.H:142
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
errorManip< error > abort(error &err)
Definition: errorManip.H:131
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.
Definition: Cloud.H:52
Foam::polyBoundaryMesh.
Output inter-processor communications stream operating on external buffer.
Definition: UOPstream.H:54
static const char nl
Definition: Ostream.H:265
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.
Definition: complex.C:47
void deleteLostParticles()
Remove lost particles from cloud and delete.
Definition: Cloud.C:111
void move(TrackCloudType &cloud, typename ParticleType::trackingData &td, const scalar trackTime)
Move the particles.
Definition: Cloud.C:142
label patchi
#define WarningInFunction
Report a warning using Foam::Warning.
const dimensionedScalar c
Speed of light in a vacuum.
label n
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:74
volScalarField & p
void writePositions() const
Write positions to <cloudName>_positions.obj file.
Definition: Cloud.C:374
fileName path(UMean.rootPath()/UMean.caseName()/functionObjects::writeFile::outputPrefix/"graphs"/UMean.instance())