All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
MPPICCloud.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) 2013-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 "MPPICCloud.H"
27 #include "PackingModel.H"
28 #include "ParticleStressModel.H"
29 #include "DampingModel.H"
30 #include "IsotropyModel.H"
31 #include "TimeScaleModel.H"
32 
33 // * * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * //
34 
35 template<class CloudType>
37 {
38  packingModel_.reset
39  (
41  (
42  this->subModelProperties(),
43  *this
44  ).ptr()
45  );
46  dampingModel_.reset
47  (
49  (
50  this->subModelProperties(),
51  *this
52  ).ptr()
53  );
54  isotropyModel_.reset
55  (
57  (
58  this->subModelProperties(),
59  *this
60  ).ptr()
61  );
62 }
63 
64 
65 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
66 
67 template<class CloudType>
69 (
70  const word& cloudName,
71  const volScalarField& rho,
72  const volVectorField& U,
73  const volScalarField& mu,
74  const dimensionedVector& g,
75  bool readFields
76 )
77 :
78  CloudType(cloudName, rho, U, mu, g, false),
79  packingModel_(nullptr),
80  dampingModel_(nullptr),
81  isotropyModel_(nullptr)
82 {
83  if (this->solution().steadyState())
84  {
86  << "MPPIC modelling not available for steady state calculations"
87  << exit(FatalError);
88  }
89 
90  if (this->solution().active())
91  {
92  setModels();
93 
94  if (readFields)
95  {
97  this->deleteLostParticles();
98  }
99  }
100 }
101 
102 
103 template<class CloudType>
105 (
107  const word& name
108 )
109 :
110  CloudType(c, name),
111  packingModel_(c.packingModel_->clone()),
112  dampingModel_(c.dampingModel_->clone()),
113  isotropyModel_(c.isotropyModel_->clone())
114 {}
115 
116 
117 template<class CloudType>
119 (
120  const fvMesh& mesh,
121  const word& name,
122  const MPPICCloud<CloudType>& c
123 )
124 :
125  CloudType(mesh, name, c),
126  packingModel_(nullptr),
127  dampingModel_(nullptr),
128  isotropyModel_(nullptr)
129 {}
130 
131 
132 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
133 
134 template<class CloudType>
136 {}
137 
138 
139 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
140 
141 template<class CloudType>
143 {
144  cloudCopyPtr_.reset
145  (
146  static_cast<MPPICCloud<CloudType>*>
147  (
148  clone(this->name() + "Copy").ptr()
149  )
150  );
151 }
152 
153 
154 template<class CloudType>
156 {
157  this->cloudReset(cloudCopyPtr_());
158  cloudCopyPtr_.clear();
159 }
160 
161 
162 template<class CloudType>
164 {
165  if (this->solution().canEvolve())
166  {
167  typename parcelType::trackingData td(*this);
168 
169  this->solve(*this, td);
170  }
171 }
172 
173 
174 template<class CloudType>
175 template<class TrackCloudType>
177 (
178  TrackCloudType& cloud,
179  typename parcelType::trackingData& td
180 )
181 {
182  // Kinematic
183  // ~~~~~~~~~
184 
185  // force calculation and tracking
186  td.part() = parcelType::trackingData::tpPredictTrack;
187  CloudType::move(cloud, td, this->db().time().deltaTValue());
188 
189 
190  // Preliminary
191  // ~~~~~~~~~~~
192 
193  // switch forces off so they are not applied in corrector steps
194  this->forces().setCalcNonCoupled(false);
195  this->forces().setCalcCoupled(false);
196 
197 
198  // Damping
199  // ~~~~~~~
200 
201  if (dampingModel_->active())
202  {
203  if (this->mesh().moving())
204  {
206  << "MPPIC damping modelling does not support moving meshes."
207  << exit(FatalError);
208  }
209 
210  // update averages
211  td.updateAverages(cloud);
212 
213  // memory allocation and eulerian calculations
214  dampingModel_->cacheFields(true);
215 
216  // calculate the damping velocity corrections without moving the parcels
217  td.part() = parcelType::trackingData::tpDampingNoTrack;
218  CloudType::move(cloud, td, this->db().time().deltaTValue());
219 
220  // correct the parcel positions and velocities
221  td.part() = parcelType::trackingData::tpCorrectTrack;
222  CloudType::move(cloud, td, this->db().time().deltaTValue());
223 
224  // finalise and free memory
225  dampingModel_->cacheFields(false);
226  }
227 
228 
229  // Packing
230  // ~~~~~~~
231 
232  if (packingModel_->active())
233  {
234  if (this->mesh().moving())
235  {
237  << "MPPIC packing modelling does not support moving meshes."
238  << exit(FatalError);
239  }
240 
241  // same procedure as for damping
242  td.updateAverages(cloud);
243  packingModel_->cacheFields(true);
244  td.part() = parcelType::trackingData::tpPackingNoTrack;
245  CloudType::move(cloud, td, this->db().time().deltaTValue());
246  td.part() = parcelType::trackingData::tpCorrectTrack;
247  CloudType::move(cloud, td, this->db().time().deltaTValue());
248  packingModel_->cacheFields(false);
249  }
250 
251 
252  // Isotropy
253  // ~~~~~~~~
254 
255  if (isotropyModel_->active())
256  {
257  // update averages
258  td.updateAverages(cloud);
259 
260  // apply isotropy model
261  isotropyModel_->calculate();
262  }
263 
264 
265  // Final
266  // ~~~~~
267 
268  // update cell occupancy
269  this->updateCellOccupancy();
270 
271  // switch forces back on
272  this->forces().setCalcNonCoupled(true);
273  this->forces().setCalcCoupled(this->solution().coupled());
274 }
275 
276 
277 template<class CloudType>
279 {
280  CloudType::info();
281 
282  tmp<volScalarField> alpha = this->theta();
283 
284  const scalar alphaMin = gMin(alpha().primitiveField());
285  const scalar alphaMax = gMax(alpha().primitiveField());
286 
287  Info<< " Min cell volume fraction = " << alphaMin << endl;
288  Info<< " Max cell volume fraction = " << alphaMax << endl;
289 
290  if (alphaMax < small)
291  {
292  return;
293  }
294 
295  scalar nMin = great;
296 
297  forAll(this->mesh().cells(), celli)
298  {
299  const label n = this->cellOccupancy()[celli].size();
300 
301  if (n > 0)
302  {
303  const scalar nPack = n*alphaMax/alpha()[celli];
304 
305  if (nPack < nMin)
306  {
307  nMin = nPack;
308  }
309  }
310  }
311 
312  reduce(nMin, minOp<scalar>());
313 
314  Info<< " Min dense number of parcels = " << nMin << endl;
315 }
316 
317 
318 // ************************************************************************* //
virtual ~MPPICCloud()
Destructor.
Definition: MPPICCloud.C:135
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
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
DSMCCloud< dsmcParcel > CloudType
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
void readFields(const typename GeoFieldType::Mesh &mesh, const IOobjectList &objects, const HashSet< word > &selectedFields, LIFOStack< regIOobject *> &storedObjects)
Read the selected GeometricFields of the specified type.
Definition: ReadFields.C:244
error FatalError
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Type gMin(const FieldField< Field, Type > &f)
void storeState()
Store the current cloud state.
Definition: MPPICCloud.C:142
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
Base class for packing models.
Definition: MPPICCloud.H:53
autoPtr< DampingModel< MPPICCloud< CloudType > > > dampingModel_
Damping model.
Definition: MPPICCloud.H:103
void motion(TrackCloudType &cloud, typename parcelType::trackingData &td)
Particle motion.
Definition: MPPICCloud.C:177
rhoEqn solve()
MPPICCloud(const word &cloudName, const volScalarField &rho, const volVectorField &U, const volScalarField &mu, const dimensionedVector &g, bool readFields=true)
Construct given carrier gas fields.
Definition: MPPICCloud.C:69
void info()
I-O.
Definition: MPPICCloud.C:278
T clone(const T &t)
Definition: List.H:54
dynamicFvMesh & mesh
const cellShapeList & cells
autoPtr< BasicCompressibleMomentumTransportModel > New(const volScalarField &rho, const volVectorField &U, const surfaceScalarField &phi, const typename BasicCompressibleMomentumTransportModel::transportModel &transport)
autoPtr< IsotropyModel< MPPICCloud< CloudType > > > isotropyModel_
Exchange model.
Definition: MPPICCloud.H:107
A class for handling words, derived from string.
Definition: word.H:59
Type gMax(const FieldField< Field, Type > &f)
Adds MPPIC modelling to kinematic clouds.
Definition: MPPICCloud.H:66
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
Base class for collisional return-to-isotropy models.
Definition: MPPICCloud.H:59
dimensionedScalar alphaMax(laminarTransport.lookup("alphaMax"))
const List< DynamicList< molecule * > > & cellOccupancy
void setModels()
Set cloud sub-models.
Definition: MPPICCloud.C:36
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
messageStream Info
label n
Selector class for relaxation factors, solver type and solution.
Definition: solution.H:48
void evolve()
Evolve the cloud.
Definition: MPPICCloud.C:163
volScalarField alpha(IOobject("alpha", runTime.timeName(), mesh, IOobject::READ_IF_PRESENT, IOobject::AUTO_WRITE), lambda *max(Ua &U, zeroSensitivity))
A class for managing temporary objects.
Definition: PtrList.H:53
Base class for collisional damping models.
Definition: MPPICCloud.H:56
autoPtr< PackingModel< MPPICCloud< CloudType > > > packingModel_
Packing model.
Definition: MPPICCloud.H:99
void restoreState()
Reset the current cloud to the previously stored state.
Definition: MPPICCloud.C:155