RecycleInteraction.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2020-2022 OpenCFD Ltd.
9 -------------------------------------------------------------------------------
10 License
11  This file is part of OpenFOAM.
12 
13  OpenFOAM is free software: you can redistribute it and/or modify it
14  under the terms of the GNU General Public License as published by
15  the Free Software Foundation, either version 3 of the License, or
16  (at your option) any later version.
17 
18  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
20  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
21  for more details.
22 
23  You should have received a copy of the GNU General Public License
24  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
25 
26 \*---------------------------------------------------------------------------*/
27 
28 #include "RecycleInteraction.H"
29 
30 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
31 
32 template<class CloudType>
34 {
36 
37  forAll(nRemoved_, i)
38  {
39  const word& outPatchName = recyclePatches_[i].first();
40 
41  forAll(nRemoved_[i], injectori)
42  {
43  const word suffix = Foam::name(injectori);
44  this->writeTabbed(os, outPatchName + "_nRemoved_" + suffix);
45  this->writeTabbed(os, outPatchName + "_massRemoved_" + suffix);
46  }
47 
48  const word& inPatchName = recyclePatches_[i].second();
49 
50  forAll(nInjected_[i], injectori)
51  {
52  const word suffix = Foam::name(injectori);
53  this->writeTabbed(os, inPatchName + "_nInjected_" + suffix);
54  this->writeTabbed(os, inPatchName + "_massInjected_" + suffix);
55  }
56  }
57 }
58 
59 
60 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
61 
62 template<class CloudType>
64 (
65  const dictionary& dict,
67 )
68 :
70  mesh_(cloud.mesh()),
71  recyclePatches_(this->coeffDict().lookup("recyclePatches")),
72  recyclePatchesIds_(recyclePatches_.size()),
73  recycledParcels_(recyclePatches_.size()),
74  nRemoved_(recyclePatches_.size()), // per patch the no. of parcels
75  massRemoved_(nRemoved_.size()),
76  nInjected_(nRemoved_.size()),
77  massInjected_(nRemoved_.size()),
78  injectionPatchPtr_(nRemoved_.size()),
79  recycleFraction_
80  (
81  this->coeffDict().template getCheck<scalar>
82  (
83  "recycleFraction",
85  )
86  ),
87  outputByInjectorId_
88  (
89  this->coeffDict().getOrDefault("outputByInjectorId", false)
90  )
91 {
92  // Determine the number of injectors and the injector mapping
93  label nInjectors = 0;
95  {
96  for (const auto& inj : cloud.injectors())
97  {
98  injIdToIndex_.insert(inj.injectorID(), ++nInjectors);
99  }
100  }
101 
102  // The normal case, and safety if injector mapping was somehow null
103  if (injIdToIndex_.empty())
104  {
105  nInjectors = 1;
106  }
107 
108  forAll(nRemoved_, i)
109  {
110  // Create injection helper for each inflow patch
112  (
113  i,
114  new patchInjectionBase(mesh_, recyclePatches_[i].second())
115  );
116 
117  // Mappings from patch names to patch IDs
118  recyclePatchesIds_[i].first() =
120  recyclePatchesIds_[i].second() =
122 
123  // Storage for reporting
124  nRemoved_[i].setSize(nInjectors, Zero);
125  massRemoved_[i].setSize(nInjectors, Zero);
126  nInjected_[i].setSize(nInjectors, Zero);
127  massInjected_[i].setSize(nInjectors, Zero);
128  }
129 }
130 
131 
132 template<class CloudType>
134 (
135  const RecycleInteraction<CloudType>& pim
136 )
137 :
138  PatchInteractionModel<CloudType>(pim),
139  mesh_(pim.mesh_),
140  recyclePatches_(pim.recyclePatches_),
141  recyclePatchesIds_(pim.recyclePatchesIds_),
142  nRemoved_(pim.nRemoved_),
143  massRemoved_(pim.massRemoved_),
144  nInjected_(pim.nInjected_),
145  massInjected_(pim.massInjected_),
146  injIdToIndex_(pim.injIdToIndex_),
147  injectionPatchPtr_(),
148  recycleFraction_(pim.recycleFraction_),
149  outputByInjectorId_(pim.outputByInjectorId_)
150 {}
151 
152 
153 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
154 
155 template<class CloudType>
157 (
158  typename CloudType::parcelType& p,
159  const polyPatch& pp,
160  bool& keepParticle
161 )
162 {
163  // Injector ID
164  const label idx =
165  (
166  injIdToIndex_.size()
167  ? injIdToIndex_.lookup(p.typeId(), 0)
168  : 0
169  );
170 
171  // See if this patch is designated an outflow patch
172  label addri = -1;
173  forAll(recyclePatchesIds_, i)
174  {
175  if (recyclePatchesIds_[i].first() == pp.index())
176  {
177  addri = i;
178  break;
179  }
180  }
181 
182  if (addri == -1)
183  {
184  // Not processing this outflow patch
185  keepParticle = true;
186  return false;
187  }
188 
189  // Flag to remove current parcel and copy to local storage
190  keepParticle = false;
191  recycledParcels_[addri].append
192  (
193  static_cast<parcelType*>(p.clone().ptr())
194  );
195 
196  ++nRemoved_[addri][idx];
197  massRemoved_[addri][idx] += p.nParticle()*p.mass();
198 
199  return true;
200 }
201 
202 
203 template<class CloudType>
205 {
206  if (Pstream::parRun())
207  {
208  // See comments in Cloud::move() about transfer particles handling
209 
210  // Allocate transfer buffers
212 
213  // Cache of opened UOPstream wrappers
214  PtrList<UOPstream> UOPstreamPtrs(Pstream::nProcs());
215 
216  auto& rnd = this->owner().rndGen();
217 
218  forAll(recycledParcels_, addri)
219  {
220  auto& patchParcels = recycledParcels_[addri];
221  auto& injectionPatch = injectionPatchPtr_[addri];
222 
223  for (parcelType& p : patchParcels)
224  {
225  // Choose a random location to insert the parcel
226  const scalar fraction01 = rnd.template sample01<scalar>();
227 
228  // Identify the processor that owns the location
229  const label toProci = injectionPatch.whichProc(fraction01);
230 
231  // Get/create output stream
232  auto* osptr = UOPstreamPtrs.get(toProci);
233  if (!osptr)
234  {
235  osptr = new UOPstream(toProci, pBufs);
236  UOPstreamPtrs.set(toProci, osptr);
237  }
238 
239  // Tuple: (address fraction particle)
240  (*osptr) << addri << fraction01 << p;
241 
242  // Can now remove from list and delete
243  delete(patchParcels.remove(&p));
244  }
245  }
246 
247  pBufs.finishedSends();
248 
249  // Not looping, so no early exit needed
250  //
251  // if (!returnReduceOr(pBufs.hasRecvData()))
252  // {
253  // // No parcels to recycle
254  // return;
255  // }
256 
257  // Retrieve from receive buffers
258  for (const int proci : pBufs.allProcs())
259  {
260  if (pBufs.recvDataCount(proci))
261  {
262  UIPstream is(proci, pBufs);
263 
264  // Read out each (address fraction particle) tuple
265  while (!is.eof())
266  {
267  const label addri = pTraits<label>(is);
268  const scalar fraction01 = pTraits<scalar>(is);
269  auto* newp = new parcelType(this->owner().mesh(), is);
270 
271  // Parcel to be recycled
272  vector newPosition;
273  label cellOwner;
274  label dummy;
275  injectionPatchPtr_[addri].setPositionAndCell
276  (
277  mesh_,
278  fraction01,
279  this->owner().rndGen(),
280  newPosition,
281  cellOwner,
282  dummy,
283  dummy
284  );
285  newp->relocate(newPosition, cellOwner);
286  newp->nParticle() *= recycleFraction_;
287 
288  // Assume parcel velocity is same as the carrier velocity
289  newp->U() = this->owner().U()[cellOwner];
290 
291  // Injector ID
292  const label idx =
293  (
294  injIdToIndex_.size()
295  ? injIdToIndex_.lookup(newp->typeId(), 0)
296  : 0
297  );
298 
299  ++nInjected_[addri][idx];
300  massInjected_[addri][idx] += newp->nParticle()*newp->mass();
301 
302  this->owner().addParticle(newp);
303  }
304  }
305  }
306  }
307  else
308  {
309  forAll(recycledParcels_, addri)
310  {
311  for (parcelType& p : recycledParcels_[addri])
312  {
313  parcelType* newp = recycledParcels_[addri].remove(&p);
314 
315  // Parcel to be recycled
316  vector newPosition;
317  label cellOwner;
318  label dummy;
319  injectionPatchPtr_[addri].setPositionAndCell
320  (
321  mesh_,
322  this->owner().rndGen(),
323  newPosition,
324  cellOwner,
325  dummy,
326  dummy
327  );
328 
329  // Update parcel properties
330  newp->relocate(newPosition, cellOwner);
331  newp->nParticle() *= recycleFraction_;
332 
333  // Assume parcel velocity is same as the carrier velocity
334  newp->U() = this->owner().U()[cellOwner];
335 
336  // Injector ID
337  const label idx =
338  (
339  injIdToIndex_.size()
340  ? injIdToIndex_.lookup(newp->typeId(), 0)
341  : 0
342  );
343  ++nInjected_[addri][idx];
344  massInjected_[addri][idx] += newp->nParticle()*newp->mass();
345 
346  this->owner().addParticle(newp);
347  }
348  }
349  }
350 }
351 
352 
353 template<class CloudType>
355 {
357 
358  labelListList npr0(nRemoved_.size());
359  scalarListList mpr0(massRemoved_.size());
360  labelListList npi0(nInjected_.size());
361  scalarListList mpi0(massInjected_.size());
362 
363  forAll(nRemoved_, patchi)
364  {
365  const label lsd = nRemoved_[patchi].size();
366  npr0[patchi].setSize(lsd, Zero);
367  mpr0[patchi].setSize(lsd, Zero);
368  npi0[patchi].setSize(lsd, Zero);
369  mpi0[patchi].setSize(lsd, Zero);
370  }
371 
372  this->getModelProperty("nRemoved", npr0);
373  this->getModelProperty("massRemoved", mpr0);
374  this->getModelProperty("nInjected", npi0);
375  this->getModelProperty("massInjected", mpi0);
376 
377  // Accumulate current data
378  labelListList npr(nRemoved_);
379 
380  forAll(npr, i)
381  {
382  Pstream::listCombineGather(npr[i], plusEqOp<label>());
383  npr[i] = npr[i] + npr0[i];
384  }
385 
386  scalarListList mpr(massRemoved_);
387  forAll(mpr, i)
388  {
389  Pstream::listCombineGather(mpr[i], plusEqOp<scalar>());
390  mpr[i] = mpr[i] + mpr0[i];
391  }
392 
393  labelListList npi(nInjected_);
394  forAll(npi, i)
395  {
396  Pstream::listCombineGather(npi[i], plusEqOp<label>());
397  npi[i] = npi[i] + npi0[i];
398  }
399 
400  scalarListList mpi(massInjected_);
401  forAll(mpi, i)
402  {
403  Pstream::listCombineGather(mpi[i], plusEqOp<scalar>());
404  mpi[i] = mpi[i] + mpi0[i];
405  }
406 
407  if (injIdToIndex_.size())
408  {
409  // Since injIdToIndex_ is a one-to-one mapping (starting as zero),
410  // can simply invert it.
411  labelList indexToInjector(injIdToIndex_.size());
412  forAllConstIters(injIdToIndex_, iter)
413  {
414  indexToInjector[iter.val()] = iter.key();
415  }
416 
417  forAll(npr, i)
418  {
419  const word& outPatchName = recyclePatches_[i].first();
420 
421  Log_<< " Parcel fate: patch " << outPatchName
422  << " (number, mass)" << nl;
423 
424  forAll(mpr[i], indexi)
425  {
426  Log_<< " - removed (injector " << indexToInjector[indexi]
427  << ") = " << npr[i][indexi]
428  << ", " << mpr[i][indexi] << nl;
429 
430  this->file()
431  << tab << npr[i][indexi] << tab << mpr[i][indexi];
432  }
433 
434  const word& inPatchName = recyclePatches_[i].second();
435 
436  Log_<< " Parcel fate: patch " << inPatchName
437  << " (number, mass)" << nl;
438 
439  forAll(mpi[i], indexi)
440  {
441  Log_<< " - injected (injector " << indexToInjector[indexi]
442  << ") = " << npi[i][indexi]
443  << ", " << mpi[i][indexi] << nl;
444  this->file()
445  << tab << npi[i][indexi] << tab << mpi[i][indexi];
446  }
447  }
448 
449  this->file() << endl;
450  }
451  else
452  {
453  forAll(npr, i)
454  {
455  const word& outPatchName = recyclePatches_[i].first();
456 
457  Log_<< " Parcel fate: patch " << outPatchName
458  << " (number, mass)" << nl
459  << " - removed = " << npr[i][0] << ", " << mpr[i][0]
460  << nl;
461 
462  this->file()
463  << tab << npr[i][0] << tab << mpr[i][0];
464  }
465 
466  forAll(npi, i)
467  {
468  const word& inPatchName = recyclePatches_[i].second();
469 
470  Log_<< " Parcel fate: patch " << inPatchName
471  << " (number, mass)" << nl
472  << " - injected = " << npi[i][0] << ", " << mpi[i][0]
473  << nl;
474 
475  this->file()
476  << tab << npi[i][0] << tab << mpi[i][0];
477  }
478 
479  this->file() << endl;
480  }
481 
482  if (this->writeTime())
483  {
484  this->setModelProperty("nRemoved", npr);
485  this->setModelProperty("massRemoved", mpr);
486  this->setModelProperty("nInjected", npi);
487  this->setModelProperty("massInjected", mpi);
488 
489  nRemoved_ = Zero;
490  massRemoved_ = Zero;
491  nInjected_ = Zero;
492  massInjected_ = Zero;
493  }
494 }
495 
496 
497 // ************************************************************************* //
label findPatchID(const word &patchName, const bool allowNotFound=true) const
Find patch index given a name, return -1 if not found.
dictionary dict
DSMCCloud< dsmcParcel > CloudType
List< List< scalar > > massInjected_
Mass of parcels injected.
A list of keyword definitions, which are a keyword followed by a number of values (eg...
Definition: dictionary.H:129
constexpr char nl
The newline &#39;\n&#39; character (0x0a)
Definition: Ostream.H:50
Random rndGen
Definition: createFields.H:23
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:531
PtrList< patchInjectionBase > injectionPatchPtr_
Injection patch pointer.
A traits class, which is primarily used for primitives and vector-space.
Definition: pTraits.H:75
static bool & parRun() noexcept
Test if this a parallel run.
Definition: UPstream.H:1049
List< List< scalar > > massRemoved_
Mass of parcels removed.
List< List< label > > nRemoved_
Number of parcels removed.
constexpr char tab
The tab &#39;\t&#39; character(0x09)
Definition: Ostream.H:49
const T * get(const label i) const
Return const pointer to element (can be nullptr), or nullptr for out-of-range access (ie...
Definition: UPtrListI.H:134
Lookup type of boundary radiation properties.
Definition: lookup.H:57
List< labelList > labelListList
List of labelList.
Definition: labelList.H:38
virtual void info()
Write patch interaction info.
Templated patch interaction model class.
bool outputByInjectorId_
Flag to output escaped/mass particles sorted by injectorID.
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:421
bool insert(const Key &key, const T &obj)
Copy insert a new entry, not overwriting existing entries.
Definition: HashTableI.H:152
virtual void postEvolve()
Post-evolve hook.
Input inter-processor communications stream using MPI send/recv etc. - operating on external buffer...
Definition: UIPstream.H:287
const fvMesh & mesh_
Reference to mesh.
List< Pair< word > > recyclePatches_
Outlet-inlet patch pair to apply parcel recycling.
static label nProcs(const label communicator=worldComm)
Number of ranks in parallel run (for given communicator). It is 1 for serial run. ...
Definition: UPstream.H:1065
Map< label > injIdToIndex_
Injector ID to local index map.
void setSize(const label n)
Alias for resize()
Definition: List.H:316
dynamicFvMesh & mesh
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for expressions::valueTypeCode::INVALID.
Definition: exprTraits.C:127
#define Log_
Report write to Foam::Info if the class log switch is true.
List< scalarList > scalarListList
List of scalarList.
Definition: scalarList.H:35
const polyBoundaryMesh & boundaryMesh() const noexcept
Return boundary mesh.
Definition: polyMesh.H:608
A class for handling words, derived from Foam::string.
Definition: word.H:63
A cloud is a registry collection of lagrangian particles.
Definition: cloud.H:53
void finishedSends(const bool wait=true)
Mark the send phase as being finished.
Vector< scalar > vector
Definition: vector.H:57
const T * set(const label i) const
Return const pointer to element (can be nullptr), or nullptr for out-of-range access (ie...
Definition: PtrList.H:159
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:290
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:56
Output inter-processor communications stream using MPI send/recv etc. - operating on external buffer...
Definition: UOPstream.H:395
bool empty() const noexcept
True if the hash table is empty.
Definition: HashTable.H:337
A Vector of values with scalar precision, where scalar is float/double depending on the compilation f...
label recvDataCount(const label proci) const
Number of unconsumed receive bytes for the specified processor. Must call finishedSends() or other fi...
OBJstream os(runTime.globalPath()/outputName)
Buffers for inter-processor communications streams (UOPstream, UIPstream).
Represents 0/1 range or concept. Used for tagged dispatch or clamping.
Definition: pTraits.H:65
virtual void writeFileHeader(Ostream &os)
Output file header information.
tmp< GeometricField< Type, PatchField, GeoMesh > > clone() const
Clone.
A list of pointers to objects of type <T>, with allocation/deallocation management of the pointers...
Definition: List.H:55
virtual void info()
Write patch interaction info.
List< List< label > > nInjected_
Number of parcels injected.
virtual bool correct(typename CloudType::parcelType &p, const polyPatch &pp, bool &keepParticle)
Apply velocity correction.
"nonBlocking" : (MPI_Isend, MPI_Irecv)
RecycleInteraction(const dictionary &dict, CloudType &cloud)
Construct from dictionary.
static void listCombineGather(const List< commsStruct > &comms, List< T > &values, const CombineOp &cop, const int tag, const label comm)
List< label > labelList
A List of labels.
Definition: List.H:62
volScalarField & p
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:69
List< Pair< label > > recyclePatchesIds_
Patch IDs of recyclePatches.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:67
bool eof() const noexcept
True if end of input seen.
Definition: IOstream.H:289
UPstream::rangeType allProcs() const noexcept
Range of ranks indices associated with PstreamBuffers.
uindirectPrimitivePatch pp(UIndirectList< face >(mesh.faces(), faceLabels), mesh.points())
forAllConstIters(mixture.phases(), phase)
Definition: pEqn.H:28
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:127