MPPICCloud.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 | Copyright (C) 2013-2016 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_(NULL),
80  dampingModel_(NULL),
81  isotropyModel_(NULL)
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  }
98  }
99 }
100 
101 
102 template<class CloudType>
104 (
106  const word& name
107 )
108 :
109  CloudType(c, name),
110  packingModel_(c.packingModel_->clone()),
111  dampingModel_(c.dampingModel_->clone()),
112  isotropyModel_(c.isotropyModel_->clone())
113 {}
114 
115 
116 template<class CloudType>
118 (
119  const fvMesh& mesh,
120  const word& name,
121  const MPPICCloud<CloudType>& c
122 )
123 :
124  CloudType(mesh, name, c),
125  packingModel_(NULL),
126  dampingModel_(NULL),
127  isotropyModel_(NULL)
128 {}
129 
130 
131 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
132 
133 template<class CloudType>
135 {}
136 
137 
138 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
139 
140 template<class CloudType>
142 {
143  cloudCopyPtr_.reset
144  (
145  static_cast<MPPICCloud<CloudType>*>
146  (
147  clone(this->name() + "Copy").ptr()
148  )
149  );
150 }
151 
152 
153 template<class CloudType>
155 {
156  this->cloudReset(cloudCopyPtr_());
157  cloudCopyPtr_.clear();
158 }
159 
160 
161 template<class CloudType>
163 {
164  if (this->solution().canEvolve())
165  {
166  typename parcelType::template
167  TrackingData<MPPICCloud<CloudType>> td(*this);
168 
169  this->solve(td);
170  }
171 }
172 
173 
174 template<class CloudType>
175 template<class TrackData>
177 {
178  // Kinematic
179  // ~~~~~~~~~
180 
181  // force calculation and tracking
182  td.part() = TrackData::tpLinearTrack;
183  CloudType::move(td, this->db().time().deltaTValue());
184 
185 
186  // Preliminary
187  // ~~~~~~~~~~~
188 
189  // switch forces off so they are not applied in corrector steps
190  this->forces().setCalcNonCoupled(false);
191  this->forces().setCalcCoupled(false);
192 
193 
194  // Damping
195  // ~~~~~~~
196 
197  if (dampingModel_->active())
198  {
199  // update averages
200  td.updateAverages(*this);
201 
202  // memory allocation and eulerian calculations
203  dampingModel_->cacheFields(true);
204 
205  // calculate the damping velocity corrections without moving the parcels
206  td.part() = TrackData::tpDampingNoTrack;
207  CloudType::move(td, this->db().time().deltaTValue());
208 
209  // correct the parcel positions and velocities
210  td.part() = TrackData::tpCorrectTrack;
211  CloudType::move(td, this->db().time().deltaTValue());
212 
213  // finalise and free memory
214  dampingModel_->cacheFields(false);
215  }
216 
217 
218  // Packing
219  // ~~~~~~~
220 
221  if (packingModel_->active())
222  {
223  // same procedure as for damping
224  td.updateAverages(*this);
225  packingModel_->cacheFields(true);
226  td.part() = TrackData::tpPackingNoTrack;
227  CloudType::move(td, this->db().time().deltaTValue());
228  td.part() = TrackData::tpCorrectTrack;
229  CloudType::move(td, this->db().time().deltaTValue());
230  packingModel_->cacheFields(false);
231  }
232 
233 
234  // Isotropy
235  // ~~~~~~~~
236 
237  if (isotropyModel_->active())
238  {
239  // update averages
240  td.updateAverages(*this);
241 
242  // apply isotropy model
243  isotropyModel_->calculate();
244  }
245 
246 
247  // Final
248  // ~~~~~
249 
250  // update cell occupancy
251  this->updateCellOccupancy();
252 
253  // switch forces back on
254  this->forces().setCalcNonCoupled(true);
255  this->forces().setCalcCoupled(this->solution().coupled());
256 }
257 
258 
259 template<class CloudType>
261 {
262  CloudType::info();
263 
264  tmp<volScalarField> alpha = this->theta();
265 
266  const scalar alphaMin = gMin(alpha().primitiveField());
267  const scalar alphaMax = gMax(alpha().primitiveField());
268 
269  Info<< " Min cell volume fraction = " << alphaMin << endl;
270  Info<< " Max cell volume fraction = " << alphaMax << endl;
271 
272  if (alphaMax < SMALL)
273  {
274  return;
275  }
276 
277  scalar nMin = GREAT;
278 
279  forAll(this->mesh().cells(), celli)
280  {
281  const label n = this->cellOccupancy()[celli].size();
282 
283  if (n > 0)
284  {
285  const scalar nPack = n*alphaMax/alpha()[celli];
286 
287  if (nPack < nMin)
288  {
289  nMin = nPack;
290  }
291  }
292  }
293 
294  reduce(nMin, minOp<scalar>());
295 
296  Info<< " Min dense number of parcels = " << nMin << endl;
297 }
298 
299 
300 // ************************************************************************* //
virtual ~MPPICCloud()
Destructor.
Definition: MPPICCloud.C:134
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
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
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:141
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:253
Base class for packing models.
Definition: MPPICCloud.H:53
autoPtr< DampingModel< MPPICCloud< CloudType > > > dampingModel_
Damping model.
Definition: MPPICCloud.H:112
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
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
void info()
I-O.
Definition: MPPICCloud.C:260
dynamicFvMesh & mesh
const cellShapeList & cells
autoPtr< IsotropyModel< MPPICCloud< CloudType > > > isotropyModel_
Exchange model.
Definition: MPPICCloud.H:116
A class for handling words, derived from string.
Definition: word.H:59
rhoEqn solve()
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
void motion(TrackData &td)
Particle motion.
Definition: MPPICCloud.C:176
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:162
A class for managing temporary objects.
Definition: PtrList.H:54
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
Base class for collisional damping models.
Definition: MPPICCloud.H:56
autoPtr< PackingModel< MPPICCloud< CloudType > > > packingModel_
Packing model.
Definition: MPPICCloud.H:108
void restoreState()
Reset the current cloud to the previously stored state.
Definition: MPPICCloud.C:154