MPPICParcelTrackingDataI.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) 2013-2022 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 "AveragingMethod.H"
27 
28 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
29 
30 template<class ParcelType>
31 template<class TrackCloudType>
33 (
34  const TrackCloudType& cloud
35 )
36 :
37  ParcelType::trackingData(cloud),
38  volumeAverage_
39  (
40  AveragingMethod<scalar>::New
41  (
42  IOobject
43  (
44  cloud.name() + ":volumeAverage",
45  cloud.db().time().name(),
46  cloud.mesh()
47  ),
48  cloud.solution().dict(),
49  cloud.mesh()
50  )
51  ),
52  radiusAverage_
53  (
54  AveragingMethod<scalar>::New
55  (
56  IOobject
57  (
58  cloud.name() + ":radiusAverage",
59  cloud.db().time().name(),
60  cloud.mesh()
61  ),
62  cloud.solution().dict(),
63  cloud.mesh()
64  )
65  ),
66  rhoAverage_
67  (
68  AveragingMethod<scalar>::New
69  (
70  IOobject
71  (
72  cloud.name() + ":rhoAverage",
73  cloud.db().time().name(),
74  cloud.mesh()
75  ),
76  cloud.solution().dict(),
77  cloud.mesh()
78  )
79  ),
80  uAverage_
81  (
83  (
84  IOobject
85  (
86  cloud.name() + ":uAverage",
87  cloud.db().time().name(),
88  cloud.mesh()
89  ),
90  cloud.solution().dict(),
91  cloud.mesh()
92  )
93  ),
94  uSqrAverage_
95  (
96  AveragingMethod<scalar>::New
97  (
98  IOobject
99  (
100  cloud.name() + ":uSqrAverage",
101  cloud.db().time().name(),
102  cloud.mesh()
103  ),
104  cloud.solution().dict(),
105  cloud.mesh()
106  )
107  ),
108  frequencyAverage_
109  (
110  AveragingMethod<scalar>::New
111  (
112  IOobject
113  (
114  cloud.name() + ":frequencyAverage",
115  cloud.db().time().name(),
116  cloud.mesh()
117  ),
118  cloud.solution().dict(),
119  cloud.mesh()
120  )
121  ),
122  massAverage_
123  (
124  AveragingMethod<scalar>::New
125  (
126  IOobject
127  (
128  cloud.name() + ":massAverage",
129  cloud.db().time().name(),
130  cloud.mesh()
131  ),
132  cloud.solution().dict(),
133  cloud.mesh()
134  )
135  )
136 {}
137 
138 
139 template<class ParcelType>
140 template<class TrackCloudType>
142 (
143  const TrackCloudType& cloud
144 )
145 {
146  // zero the sums
147  volumeAverage_() = 0;
148  radiusAverage_() = 0;
149  rhoAverage_() = 0;
150  uAverage_() = Zero;
151  uSqrAverage_() = 0;
152  frequencyAverage_() = 0;
153  massAverage_() = 0;
154 
155  // temporary weights
156  autoPtr<AveragingMethod<scalar>> weightAveragePtr
157  (
159  (
160  IOobject
161  (
162  cloud.name() + ":weightAverage",
163  cloud.db().time().name(),
164  cloud.mesh()
165  ),
166  cloud.solution().dict(),
167  cloud.mesh()
168  )
169  );
170  AveragingMethod<scalar>& weightAverage = weightAveragePtr();
171 
172  // averaging sums
173  forAllConstIter(typename TrackCloudType, cloud, iter)
174  {
175  const typename TrackCloudType::parcelType& p = iter();
176  const tetIndices tetIs = p.currentTetIndices(cloud.mesh());
177 
178  const scalar m = p.nParticle()*p.mass();
179 
180  volumeAverage_->add(p.coordinates(), tetIs, p.nParticle()*p.volume());
181  rhoAverage_->add(p.coordinates(), tetIs, m*p.rho());
182  uAverage_->add(p.coordinates(), tetIs, m*p.U());
183  massAverage_->add(p.coordinates(), tetIs, m);
184  }
185  volumeAverage_->average();
186  massAverage_->average();
187  rhoAverage_->average(massAverage_);
188  uAverage_->average(massAverage_);
189 
190  // squared velocity deviation
191  forAllConstIter(typename TrackCloudType, cloud, iter)
192  {
193  const typename TrackCloudType::parcelType& p = iter();
194  const tetIndices tetIs = p.currentTetIndices(cloud.mesh());
195 
196  const vector u = uAverage_->interpolate(p.coordinates(), tetIs);
197 
198  uSqrAverage_->add
199  (
200  p.coordinates(),
201  tetIs,
202  p.nParticle()*p.mass()*magSqr(p.U() - u)
203  );
204  }
205  uSqrAverage_->average(massAverage_);
206 
207  // sauter mean radius
208  radiusAverage_() = volumeAverage_();
209  weightAverage = 0;
210  forAllConstIter(typename TrackCloudType, cloud, iter)
211  {
212  const typename TrackCloudType::parcelType& p = iter();
213  const tetIndices tetIs = p.currentTetIndices(cloud.mesh());
214 
215  weightAverage.add
216  (
217  p.coordinates(),
218  tetIs,
219  p.nParticle()*pow(p.volume(), 2.0/3.0)
220  );
221  }
222  weightAverage.average();
223  radiusAverage_->average(weightAverage);
224 
225  // collision frequency
226  weightAverage = 0;
227  forAllConstIter(typename TrackCloudType, cloud, iter)
228  {
229  const typename TrackCloudType::parcelType& p = iter();
230  const tetIndices tetIs = p.currentTetIndices(cloud.mesh());
231 
232  const scalar a = volumeAverage_->interpolate(p.coordinates(), tetIs);
233  const scalar r = radiusAverage_->interpolate(p.coordinates(), tetIs);
234  const vector u = uAverage_->interpolate(p.coordinates(), tetIs);
235 
236  const scalar f = 0.75*a/pow3(r)*sqr(0.5*p.d() + r)*mag(p.U() - u);
237 
238  frequencyAverage_->add(p.coordinates(), tetIs, p.nParticle()*f*f);
239 
240  weightAverage.add(p.coordinates(), tetIs, p.nParticle()*f);
241  }
242  frequencyAverage_->average(weightAverage);
243 }
244 
245 
246 // ************************************************************************* //
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
Base class for lagrangian averaging methods.
virtual void add(const barycentric &coordinates, const tetIndices &tetIs, const Type &value)=0
Member Functions.
virtual void average()
Calculate the average.
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
const objectRegistry & db() const
Return the local objectRegistry.
Definition: IOobject.C:312
const word & name() const
Return name.
Definition: IOobject.H:310
trackingData(const TrackCloudType &cloud)
Constructors.
void updateAverages(const TrackCloudType &cloud)
Update the MPPIC averages.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: autoPtr.H:51
A cloud is a collection of lagrangian particles.
Definition: cloud.H:55
const word & name() const
Return const reference to name.
const Time & time() const
Return time.
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:51
Storage and named access for the indices of a tet which is part of the decomposition of a cell.
Definition: tetIndices.H:82
autoPtr< CompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const viscosity &viscosity)
static const zero Zero
Definition: zero.H:97
dimensionedSymmTensor sqr(const dimensionedVector &dv)
dimensionedScalar pow3(const dimensionedScalar &ds)
word name(const bool)
Return a word representation of a bool.
Definition: boolIO.C:39
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
dimensioned< scalar > mag(const dimensioned< Type > &)
dimensioned< scalar > magSqr(const dimensioned< Type > &)
labelList f(nPoints)
dictionary dict
volScalarField & p