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 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",
84  scalarMinMax::zero_one()
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;
94  if (outputByInjectorId_)
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
111  injectionPatchPtr_.set
112  (
113  i,
114  new patchInjectionBase(mesh_, recyclePatches_[i].second())
115  );
116 
117  // Mappings from patch names to patch IDs
118  recyclePatchesIds_[i].first() =
119  mesh_.boundaryMesh().findPatchID(recyclePatches_[i].first());
120  recyclePatchesIds_[i].second() =
121  mesh_.boundaryMesh().findPatchID(recyclePatches_[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 (
136 )
137 :
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  // Relocate the recycled parcels into slots for each receiving processor
209  List<IDLList<parcelType>> transferParcels(Pstream::nProcs());
210  List<DynamicList<scalar>> fractions(Pstream::nProcs());
211  List<DynamicList<label>> patchAddr(Pstream::nProcs());
212 
213  auto& rnd = this->owner().rndGen();
214 
215  forAll(recycledParcels_, addri)
216  {
217  auto& patchParcels = recycledParcels_[addri];
218  auto& injectionPatch = injectionPatchPtr_[addri];
219 
220  forAllIters(patchParcels, pIter)
221  {
222  // Choose a random location to insert the parcel
223  const scalar fraction01 = rnd.template sample01<scalar>();
224 
225  // Identify the processor that owns the location
226  const label toProci = injectionPatch.whichProc(fraction01);
227 
228  // Store info in slot for target processor
229  transferParcels[toProci].append(patchParcels.remove(pIter));
230  fractions[toProci].append(fraction01);
231  patchAddr[toProci].append(addri);
232  }
233  }
234 
235  // Set-up the sends
236  PstreamBuffers pBufs(Pstream::commsTypes::nonBlocking);
237 
238  // Clear transfer buffers
239  pBufs.clear();
240 
241  // Stream into send buffers
242  forAll(transferParcels, proci)
243  {
244  if (transferParcels[proci].size())
245  {
246  UOPstream particleStream(proci, pBufs);
247 
248  particleStream
249  << transferParcels[proci]
250  << fractions[proci]
251  << patchAddr[proci];
252 
253  transferParcels[proci].clear();
254  }
255  }
256 
257  // Start sending. Sets number of bytes transferred
258  labelList allNTrans(Pstream::nProcs());
259  pBufs.finishedSends(allNTrans);
260  bool transferred = false;
261  for (const label n : allNTrans)
262  {
263  if (n)
264  {
265  transferred = true;
266  break;
267  }
268  }
269  reduce(transferred, orOp<bool>());
270  if (!transferred)
271  {
272  // No parcels to transfer
273  return;
274  }
275 
276  // Retrieve from receive buffers
277  for (label proci = 0; proci < Pstream::nProcs(); ++proci)
278  {
279  if (allNTrans[proci])
280  {
281  UIPstream particleStream(proci, pBufs);
282  IDLList<parcelType> newParticles
283  (
284  particleStream,
285  typename parcelType::iNew(this->owner().mesh())
286  );
287  scalarList fractions(particleStream);
288  labelList patchAddr(particleStream);
289 
290  label parceli = 0;
291  for (parcelType& newp : newParticles)
292  {
293  // Parcel to be recycled
294  vector newPosition;
295  label cellOwner;
296  label dummy;
297  const label addri = patchAddr[parceli];
298  injectionPatchPtr_[addri].setPositionAndCell
299  (
300  mesh_,
301  fractions[parceli],
302  this->owner().rndGen(),
303  newPosition,
304  cellOwner,
305  dummy,
306  dummy
307  );
308  newp.relocate(newPosition, cellOwner);
309  newp.U() = this->owner().U()[cellOwner];
310  newp.nParticle() *= recycleFraction_;
311 
312  const label idx =
313  (
314  injIdToIndex_.size()
315  ? injIdToIndex_.lookup(newp.typeId(), 0)
316  : 0
317  );
318  ++nInjected_[addri][idx];
319  massInjected_[addri][idx] += newp.nParticle()*newp.mass();
320 
321  this->owner().addParticle(newParticles.remove(&newp));
322  ++parceli;
323  }
324  }
325  }
326  }
327  else
328  {
329  forAll(recycledParcels_, addri)
330  {
331  forAllIters(recycledParcels_[addri], iter)
332  {
333  // Parcel to be recycled
334  vector newPosition;
335  label cellOwner;
336  label dummy;
337  injectionPatchPtr_[addri].setPositionAndCell
338  (
339  mesh_,
340  this->owner().rndGen(),
341  newPosition,
342  cellOwner,
343  dummy,
344  dummy
345  );
346 
347  // Update parcel properties
348  parcelType* newp = recycledParcels_[addri].remove(iter);
349  newp->relocate(newPosition, cellOwner);
350  newp->nParticle() *= recycleFraction_;
351 
352  // Assume parcel velocity is same as the carrier velocity
353  newp->U() = this->owner().U()[cellOwner];
354 
355  const label idx =
356  (
357  injIdToIndex_.size()
358  ? injIdToIndex_.lookup(newp->typeId(), 0)
359  : 0
360  );
361  ++nInjected_[addri][idx];
362  massInjected_[addri][idx] += newp->nParticle()*newp->mass();
363 
364  this->owner().addParticle(newp);
365  }
366  }
367  }
368 }
369 
370 
371 template<class CloudType>
373 {
375 
376  labelListList npr0(nRemoved_.size());
377  scalarListList mpr0(massRemoved_.size());
378  labelListList npi0(nInjected_.size());
379  scalarListList mpi0(massInjected_.size());
380 
381  forAll(nRemoved_, patchi)
382  {
383  label lsd = nRemoved_[patchi].size();
384  npr0[patchi].setSize(lsd, Zero);
385  mpr0[patchi].setSize(lsd, Zero);
386  npi0[patchi].setSize(lsd, Zero);
387  mpi0[patchi].setSize(lsd, Zero);
388  }
389 
390  this->getModelProperty("nRemoved", npr0);
391  this->getModelProperty("massRemoved", mpr0);
392  this->getModelProperty("nInjected", npi0);
393  this->getModelProperty("massInjected", mpi0);
394 
395  // Accumulate current data
396  labelListList npr(nRemoved_);
397 
398  forAll(npr, i)
399  {
400  Pstream::listCombineGather(npr[i], plusEqOp<label>());
401  npr[i] = npr[i] + npr0[i];
402  }
403 
404  scalarListList mpr(massRemoved_);
405  forAll(mpr, i)
406  {
407  Pstream::listCombineGather(mpr[i], plusEqOp<scalar>());
408  mpr[i] = mpr[i] + mpr0[i];
409  }
410 
411  labelListList npi(nInjected_);
412  forAll(npi, i)
413  {
414  Pstream::listCombineGather(npi[i], plusEqOp<label>());
415  npi[i] = npi[i] + npi0[i];
416  }
417 
418  scalarListList mpi(massInjected_);
419  forAll(mpi, i)
420  {
421  Pstream::listCombineGather(mpi[i], plusEqOp<scalar>());
422  mpi[i] = mpi[i] + mpi0[i];
423  }
424 
425  if (injIdToIndex_.size())
426  {
427  // Since injIdToIndex_ is a one-to-one mapping (starting as zero),
428  // can simply invert it.
429  labelList indexToInjector(injIdToIndex_.size());
430  forAllConstIters(injIdToIndex_, iter)
431  {
432  indexToInjector[iter.val()] = iter.key();
433  }
434 
435  forAll(npr, i)
436  {
437  const word& outPatchName = recyclePatches_[i].first();
438 
439  os << " Parcel fate: patch " << outPatchName
440  << " (number, mass)" << nl;
441 
442  forAll(mpr[i], indexi)
443  {
444  os << " - removed (injector " << indexToInjector[indexi]
445  << ") = " << npr[i][indexi]
446  << ", " << mpr[i][indexi] << nl;
447 
448  this->file()
449  << tab << npr[i][indexi] << tab << mpr[i][indexi];
450  }
451 
452  const word& inPatchName = recyclePatches_[i].second();
453 
454  os << " Parcel fate: patch " << inPatchName
455  << " (number, mass)" << nl;
456 
457  forAll(mpi[i], indexi)
458  {
459  os << " - injected (injector " << indexToInjector[indexi]
460  << ") = " << npi[i][indexi]
461  << ", " << mpi[i][indexi] << nl;
462  this->file()
463  << tab << npi[i][indexi] << tab << mpi[i][indexi];
464  }
465  }
466 
467  this->file() << endl;
468  }
469  else
470  {
471  forAll(npr, i)
472  {
473  const word& outPatchName = recyclePatches_[i].first();
474 
475  os << " Parcel fate: patch " << outPatchName
476  << " (number, mass)" << nl
477  << " - removed = " << npr[i][0] << ", " << mpr[i][0]
478  << nl;
479 
480  this->file()
481  << tab << npr[i][0] << tab << mpr[i][0];
482  }
483 
484  forAll(npi, i)
485  {
486  const word& inPatchName = recyclePatches_[i].second();
487 
488  os << " Parcel fate: patch " << inPatchName
489  << " (number, mass)" << nl
490  << " - injected = " << npi[i][0] << ", " << mpi[i][0]
491  << nl;
492 
493  this->file()
494  << tab << npi[i][0] << tab << mpi[i][0];
495  }
496 
497  this->file() << endl;
498  }
499 
500  if (this->writeTime())
501  {
502  this->setModelProperty("nRemoved", npr);
503  this->setModelProperty("massRemoved", mpr);
504  this->setModelProperty("nInjected", npi);
505  this->setModelProperty("massInjected", mpi);
506 
507  nRemoved_ = Zero;
508  massRemoved_ = Zero;
509  nInjected_ = Zero;
510  massInjected_ = Zero;
511  }
512 }
513 
514 
515 // ************************************************************************* //
Foam::RecycleInteraction::recyclePatchesIds_
List< Pair< label > > recyclePatchesIds_
Patch IDs of recyclePatches.
Definition: RecycleInteraction.H:162
Foam::RecycleInteraction::parcelType
CloudType::parcelType parcelType
Definition: RecycleInteraction.H:151
p
volScalarField & p
Definition: createFieldRefs.H:8
RecycleInteraction.H
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::UOPstream
Output inter-processor communications stream operating on external buffer.
Definition: UOPstream.H:57
Foam::Zero
static constexpr const zero Zero
Global zero (0)
Definition: zero.H:131
Foam::RecycleInteraction::recyclePatches_
List< Pair< word > > recyclePatches_
Outlet-inlet patch pair to apply parcel recycling.
Definition: RecycleInteraction.H:159
fractions
dictionary fractions(initialConditions.subDict("fractions"))
Foam::PstreamBuffers
Buffers for inter-processor communications streams (UOPstream, UIPstream).
Definition: PstreamBuffers.H:88
Foam::RecycleInteraction::massRemoved_
List< List< scalar > > massRemoved_
Mass of parcels removed.
Definition: RecycleInteraction.H:174
Foam::RecycleInteraction::correct
virtual bool correct(typename CloudType::parcelType &p, const polyPatch &pp, bool &keepParticle)
Apply velocity correction.
Definition: RecycleInteraction.C:157
Foam::List::append
void append(const T &val)
Append an element at the end of the list.
Definition: ListI.H:175
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::PstreamBuffers::clear
void clear()
Reset (clear) individual buffers and reset state.
Definition: PstreamBuffers.C:125
Foam::RecycleInteraction::writeFileHeader
virtual void writeFileHeader(Ostream &os)
Output file header information.
Definition: RecycleInteraction.C:33
Foam::RecycleInteraction::RecycleInteraction
RecycleInteraction(const dictionary &dict, CloudType &cloud)
Construct from dictionary.
Definition: RecycleInteraction.C:64
forAll
#define forAll(list, i)
Loop across all elements in list.
Definition: stdFoam.H:296
Foam::RecycleInteraction::mesh_
const fvMesh & mesh_
Reference to mesh.
Definition: RecycleInteraction.H:156
n
label n
Definition: TABSMDCalcMethod2.H:31
Foam::RecycleInteraction::outputByInjectorId_
bool outputByInjectorId_
Flag to output escaped/mass particles sorted by injectorID.
Definition: RecycleInteraction.H:193
Foam::PatchInteractionModel
Templated patch interaction model class.
Definition: KinematicCloud.H:89
Foam::polyPatch
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:68
Foam::RecycleInteraction::nInjected_
List< List< label > > nInjected_
Number of parcels injected.
Definition: RecycleInteraction.H:177
Foam::RecycleInteraction::postEvolve
virtual void postEvolve()
Post-evolve hook.
Definition: RecycleInteraction.C:204
Foam::DSMCCloud
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:71
Foam::PstreamBuffers::finishedSends
void finishedSends(const bool block=true)
Mark all sends as having been done.
Definition: PstreamBuffers.C:73
forAllIters
#define forAllIters(container, iter)
Iterate across all elements in the container object.
Definition: stdFoam.H:223
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::ILList
Template class for intrusive linked lists.
Definition: ILList.H:52
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
Foam::RecycleInteraction::info
virtual void info(Ostream &os)
Write patch interaction info to stream.
Definition: RecycleInteraction.C:372
reduce
reduce(hasMovingMesh, orOp< bool >())
os
OBJstream os(runTime.globalPath()/outputName)
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
Foam::patchIdentifier::index
label index() const noexcept
The index of this patch in the boundaryMesh.
Definition: patchIdentifier.H:147
Foam::RecycleInteraction
Patch interaction model to perform parcel transfer and recycle from one patch to another.
Definition: RecycleInteraction.H:145
Foam::cloud
A cloud is a registry collection of lagrangian particles.
Definition: cloud.H:57
Foam::tab
constexpr char tab
Definition: Ostream.H:403
Foam::nl
constexpr char nl
Definition: Ostream.H:404
forAllConstIters
forAllConstIters(mixture.phases(), phase)
Definition: pEqn.H:28
Foam::Vector< scalar >
Foam::patchInjectionBase
Definition: patchInjectionBase.H:64
Foam::List
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: BitOps.H:63
Foam::List::clear
void clear()
Clear the list, i.e. set size to zero.
Definition: ListI.H:116
Foam::RecycleInteraction::massInjected_
List< List< scalar > > massInjected_
Mass of parcels injected.
Definition: RecycleInteraction.H:180
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
rndGen
Random rndGen
Definition: createFields.H:23
Foam::RecycleInteraction::injIdToIndex_
Map< label > injIdToIndex_
Injector ID to local index map.
Definition: RecycleInteraction.H:184
Foam::DSMCCloud::parcelType
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:220
Foam::plusEqOp
Definition: ops.H:72
Foam::Ostream
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:56
Foam::UIPstream
Input inter-processor communications stream operating on external buffer.
Definition: UIPstream.H:56
Foam::orOp
Definition: ops.H:234
Foam::RecycleInteraction::nRemoved_
List< List< label > > nRemoved_
Number of parcels removed.
Definition: RecycleInteraction.H:171
Foam::RecycleInteraction::recycleFraction_
const scalar recycleFraction_
Parcel fraction to be recycled from outlet to inlet.
Definition: RecycleInteraction.H:190