SprayCloudI.H
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-2018 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 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
27 
28 template<class CloudType>
29 inline const Foam::SprayCloud<CloudType>&
31 {
32  return cloudCopyPtr_();
33 }
34 
35 
36 template<class CloudType>
39 {
40  return atomizationModel_;
41 }
42 
43 
44 template<class CloudType>
47 {
48  return atomizationModel_();
49 }
50 
51 
52 template<class CloudType>
55 {
56  return breakupModel_;
57 }
58 
59 
60 template<class CloudType>
63 {
64  return breakupModel_();
65 }
66 
67 
68 template<class CloudType>
70 {
71  return averageParcelMass_;
72 }
73 
74 
75 template<class CloudType>
77 (
78  const scalar fraction
79 ) const
80 {
81  if ((fraction < 0) || (fraction > 1))
82  {
84  << "fraction should be in the range 0 < fraction < 1"
85  << exit(FatalError);
86  }
87 
88  scalar distance = 0.0;
89 
90  const label nParcel = this->size();
91  globalIndex globalParcels(nParcel);
92  const label nParcelSum = globalParcels.size();
93 
94  if (nParcelSum == 0)
95  {
96  return distance;
97  }
98 
99  // lists of parcels mass and distance from initial injection point
100  List<List<scalar>> procMass(Pstream::nProcs());
101  List<List<scalar>> procDist(Pstream::nProcs());
102 
103  List<scalar>& mass = procMass[Pstream::myProcNo()];
104  List<scalar>& dist = procDist[Pstream::myProcNo()];
105 
106  mass.setSize(nParcel);
107  dist.setSize(nParcel);
108 
109  label i = 0;
110  scalar mSum = 0.0;
111  forAllConstIter(typename SprayCloud<CloudType>, *this, iter)
112  {
113  const parcelType& p = iter();
114  scalar m = p.nParticle()*p.mass();
115  scalar d = mag(p.position() - p.position0());
116  mSum += m;
117 
118  mass[i] = m;
119  dist[i] = d;
120 
121  i++;
122  }
123 
124  // calculate total mass across all processors
125  reduce(mSum, sumOp<scalar>());
126  Pstream::gatherList(procMass);
127  Pstream::gatherList(procDist);
128 
129  if (Pstream::master())
130  {
131  // flatten the mass lists
132  List<scalar> allMass(nParcelSum, 0.0);
133  SortableList<scalar> allDist(nParcelSum, 0.0);
134  for (label proci = 0; proci < Pstream::nProcs(); proci++)
135  {
137  (
138  allMass,
139  globalParcels.localSize(proci),
140  globalParcels.offset(proci)
141  ) = procMass[proci];
142 
143  // flatten the distance list
145  (
146  allDist,
147  globalParcels.localSize(proci),
148  globalParcels.offset(proci)
149  ) = procDist[proci];
150  }
151 
152  // sort allDist distances into ascending order
153  // note: allMass masses are left unsorted
154  allDist.sort();
155 
156  if (nParcelSum > 1)
157  {
158  const scalar mLimit = fraction*mSum;
159  const labelList& indices = allDist.indices();
160 
161  if (mLimit > (mSum - allMass[indices.last()]))
162  {
163  distance = allDist.last();
164  }
165  else
166  {
167  // assuming that 'fraction' is generally closer to 1 than 0,
168  // loop through in reverse distance order
169  const scalar mThreshold = (1.0 - fraction)*mSum;
170  scalar mCurrent = 0.0;
171  label i0 = 0;
172 
173  forAllReverse(indices, i)
174  {
175  label indI = indices[i];
176 
177  mCurrent += allMass[indI];
178 
179  if (mCurrent > mThreshold)
180  {
181  i0 = i;
182  break;
183  }
184  }
185 
186  if (i0 == indices.size() - 1)
187  {
188  distance = allDist.last();
189  }
190  else
191  {
192  // linearly interpolate to determine distance
193  scalar alpha = (mCurrent - mThreshold)/allMass[indices[i0]];
194  distance =
195  allDist[i0] + alpha*(allDist[i0+1] - allDist[i0]);
196  }
197  }
198  }
199  else
200  {
201  distance = allDist.first();
202  }
203  }
204 
205  Pstream::scatter(distance);
206 
207  return distance;
208 }
209 
210 
211 // ************************************************************************* //
scalar averageParcelMass() const
Return const-access to the average parcel mass.
Definition: SprayCloudI.H:69
void sort()
(stable) sort the list (if changed after construction time)
Definition: SortableList.C:104
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
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
A 1D array of objects of type <T>, where the size of the vector is known and used for subscript bound...
Definition: HashTable.H:60
A list that is sorted upon construction or when explicitly requested with the sort() method...
Definition: List.H:73
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
const AtomizationModel< SprayCloud< CloudType > > & atomization() const
Return const-access to the atomization model.
Definition: SprayCloudI.H:38
Templated atomization model class.
Definition: SprayCloud.H:47
#define forAllReverse(list, i)
Reverse loop across all elements in list.
Definition: UList.H:440
scalar penetration(const scalar fraction) const
Penetration for fraction [0-1] of the current total mass.
Definition: SprayCloudI.H:77
const SprayCloud & cloudCopy() const
Return a reference to the cloud copy.
Definition: SprayCloudI.H:30
scalar distance(const vector &p1, const vector &p2)
Definition: curveTools.C:12
T & first()
Return the first element of the list.
Definition: UListI.H:114
label offset(const label proci) const
Start of proci data.
Definition: globalIndexI.H:48
Calculates a unique integer (label so might not have enough room - 2G max) for processor + local inde...
Definition: globalIndex.H:63
A List obtained as a section of another List.
Definition: SubList.H:53
label size() const
Global sum of localSizes.
Definition: globalIndexI.H:66
forAllConstIter(PtrDictionary< phaseModel >, mixture.phases(), phase)
Definition: pEqn.H:29
const BreakupModel< SprayCloud< CloudType > > & breakup() const
Return const-access to the breakup model.
Definition: SprayCloudI.H:54
const labelList & indices() const
Return the list of sorted indices. Updated every sort.
Definition: SortableList.H:93
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:218
Templated base class for spray cloud.
Definition: SprayCloud.H:57
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
void setSize(const label)
Reset size of List.
Definition: List.C:281
dimensioned< scalar > mag(const dimensioned< Type > &)
Templated break-up model class.
Definition: SprayCloud.H:50
volScalarField & p
T & last()
Return the last element of the list.
Definition: UListI.H:128
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
label localSize() const
My local size.
Definition: globalIndexI.H:60