All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
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-2020 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(neighbourProcs.size());
178 
179  // List of destination processorPatches indices for all of the
180  // neighbour processors
181  List<DynamicList<label>> patchIndexTransferLists(neighbourProcs.size());
182 
183  // Allocate transfer buffers
184  PstreamBuffers pBufs(Pstream::commsTypes::nonBlocking);
185 
186  // Clear the global positions as there are about to change
187  globalPositionsPtr_.clear();
188 
189  // While there are particles to transfer
190  while (true)
191  {
192  // Clear the transfer lists
193  particleTransferLists = IDLList<ParticleType>();
194  forAll(patchIndexTransferLists, i)
195  {
196  patchIndexTransferLists[i].clear();
197  }
198 
199  // Loop over all particles
200  forAllIter(typename Cloud<ParticleType>, *this, pIter)
201  {
202  ParticleType& p = pIter();
203 
204  // Move the particle
205  bool keepParticle = p.move(cloud, td, trackTime);
206 
207  // If the particle is to be kept
208  // (i.e. it hasn't passed through an inlet or outlet)
209  if (keepParticle)
210  {
211  if (td.switchProcessor)
212  {
213  #ifdef FULLDEBUG
214  if
215  (
216  !Pstream::parRun()
217  || !p.onBoundaryFace()
218  || procPatchNeighbours[p.patch()] < 0
219  )
220  {
222  << "Switch processor flag is true when no parallel "
223  << "transfer is possible. This is a bug."
224  << exit(FatalError);
225  }
226  #endif
227 
228  const label patchi = p.patch();
229 
230  const label n = neighbourProcIndices
231  [
232  refCast<const processorPolyPatch>
233  (
234  pbm[patchi]
235  ).neighbProcNo()
236  ];
237 
238  p.prepareForParallelTransfer();
239 
240  particleTransferLists[n].append(this->remove(&p));
241 
242  patchIndexTransferLists[n].append
243  (
244  procPatchNeighbours[patchi]
245  );
246  }
247  }
248  else
249  {
250  deleteParticle(p);
251  }
252  }
253 
254  if (!Pstream::parRun())
255  {
256  break;
257  }
258 
259 
260  // Clear transfer buffers
261  pBufs.clear();
262 
263  // Stream into send buffers
264  forAll(particleTransferLists, i)
265  {
266  if (particleTransferLists[i].size())
267  {
268  UOPstream particleStream
269  (
270  neighbourProcs[i],
271  pBufs
272  );
273 
274  particleStream
275  << patchIndexTransferLists[i]
276  << particleTransferLists[i];
277  }
278  }
279 
280 
281  // Start sending. Sets number of bytes transferred
282  labelList allNTrans(Pstream::nProcs());
283  pBufs.finishedSends(allNTrans);
284 
285 
286  bool transferred = false;
287 
288  forAll(allNTrans, i)
289  {
290  if (allNTrans[i])
291  {
292  transferred = true;
293  break;
294  }
295  }
296  reduce(transferred, orOp<bool>());
297 
298  if (!transferred)
299  {
300  break;
301  }
302 
303  // Retrieve from receive buffers
304  forAll(neighbourProcs, i)
305  {
306  label neighbProci = neighbourProcs[i];
307 
308  label nRec = allNTrans[neighbProci];
309 
310  if (nRec)
311  {
312  UIPstream particleStream(neighbProci, pBufs);
313 
314  labelList receivePatchIndex(particleStream);
315 
316  IDLList<ParticleType> newParticles
317  (
318  particleStream,
319  typename ParticleType::iNew(polyMesh_)
320  );
321 
322  label pI = 0;
323 
324  forAllIter(typename Cloud<ParticleType>, newParticles, newpIter)
325  {
326  ParticleType& newp = newpIter();
327 
328  label patchi = procPatches[receivePatchIndex[pI++]];
329 
330  newp.correctAfterParallelTransfer(patchi, td);
331 
332  addParticle(newParticles.remove(&newp));
333  }
334  }
335  }
336  }
337 }
338 
339 
340 template<class ParticleType>
342 {
343  if (!globalPositionsPtr_.valid())
344  {
346  << "Global positions are not available. "
347  << "Cloud::storeGlobalPositions has not been called."
348  << exit(FatalError);
349  }
350 
351  // Ask for the tetBasePtIs to trigger all processors to build
352  // them, otherwise, if some processors have no particles then
353  // there is a comms mismatch.
354  polyMesh_.tetBasePtIs();
355  polyMesh_.oldCellCentres();
356 
357  const vectorField& positions = globalPositionsPtr_();
358 
359  label i = 0;
360  forAllIter(typename Cloud<ParticleType>, *this, iter)
361  {
362  iter().autoMap(positions[i], mapper);
363  ++ i;
364  }
365 }
366 
367 
368 template<class ParticleType>
370 {
371  OFstream pObj
372  (
373  this->db().time().path()/this->name() + "_positions.obj"
374  );
375 
376  forAllConstIter(typename Cloud<ParticleType>, *this, pIter)
377  {
378  const ParticleType& p = pIter();
379  pObj<< "v " << p.position().x() << " " << p.position().y() << " "
380  << p.position().z() << nl;
381  }
382 
383  pObj.flush();
384 }
385 
386 
387 template<class ParticleType>
389 {
390  // Store the global positions for later use by autoMap. It would be
391  // preferable not to need this. If the mapPolyMesh object passed to autoMap
392  // had a copy of the old mesh then the global positions could be recovered
393  // within autoMap, and this pre-processing would not be necessary.
394 
395  globalPositionsPtr_.reset(new vectorField(this->size()));
396 
397  vectorField& positions = globalPositionsPtr_();
398 
399  label i = 0;
400  forAllConstIter(typename Cloud<ParticleType>, *this, iter)
401  {
402  positions[i] = iter().position();
403  ++ i;
404  }
405 }
406 
407 
408 // * * * * * * * * * * * * * * * * IOStream operators * * * * * * * * * * * //
409 
410 #include "CloudIO.C"
411 
412 // ************************************************************************* //
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:164
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:251
volVectorField vectorField(fieldObject, mesh)
const dimensionedScalar & c
Speed of light in a vacuum.
void autoMap(const mapPolyMesh &)
Remap the cells of particles corresponding to the.
Definition: Cloud.C:341
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:388
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:125
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:178
virtual void flush()
Flush stream.
Definition: OSstream.C:207
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
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:260
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.
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:369
fileName path(UMean.rootPath()/UMean.caseName()/functionObjects::writeFile::outputPrefix/"graphs"/UMean.instance())