DSMCCloud.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-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 Class
25  Foam::DSMCCloud
26 
27 Description
28  Templated base class for dsmc cloud
29 
30 SourceFiles
31  DSMCCloudI.H
32  DSMCCloud.C
33 
34 \*---------------------------------------------------------------------------*/
35 
36 #ifndef DSMCCloud_H
37 #define DSMCCloud_H
38 
39 #include "Cloud.H"
40 #include "IOdictionary.H"
41 #include "autoPtr.H"
42 #include "Random.H"
43 #include "fvMesh.H"
44 #include "volFields.H"
45 #include "scalarIOField.H"
46 #include "barycentric.H"
47 
48 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
49 
50 namespace Foam
51 {
52 
53 // Forward declaration of classes
54 
55 template<class CloudType>
57 
58 template<class CloudType>
60 
61 template<class CloudType>
63 
64 /*---------------------------------------------------------------------------*\
65  Class DSMCCloudName Declaration
66 \*---------------------------------------------------------------------------*/
67 
69 
70 
71 /*---------------------------------------------------------------------------*\
72  Class DSMCCloud Declaration
73 \*---------------------------------------------------------------------------*/
74 
75 template<class ParcelType>
76 class DSMCCloud
77 :
78  public Cloud<ParcelType>,
79  public DSMCCloudName
80 {
81  // Private Data
82 
83  //- Cloud type - used to set the name of the parcel properties
84  // dictionary by appending "Properties"
85  const word cloudName_;
86 
87  //- References to the mesh and time databases
88  const fvMesh& mesh_;
89 
90  //- Dictionary of particle properties
91  IOdictionary particleProperties_;
92 
93  //- A list of unique instances of molecule types in the
94  // simulation. The position of an entry in the list maps to
95  // the label identifying the typeId, i.e. where typeIdList_ =
96  // (N2 O2 CO2) N2 has typeId label = 0, O2 = 1, CO2 = 2.
97  List<word> typeIdList_;
98 
99  //- Number of real atoms/molecules represented by a parcel
100  scalar nParticle_;
101 
102  //- A data structure holding which particles are in which cell
103  List<DynamicList<ParcelType*>> cellOccupancy_;
104 
105  //- A field holding the value of (sigmaT * cR)max for each
106  // cell (see Bird p220). Initialised with the parcels,
107  // updated as required, and read in on start/restart.
108  volScalarField sigmaTcRMax_;
109 
110  //- A field holding the remainder from the previous collision selections
111  volScalarField::Internal collisionSelectionRemainder_;
112 
113  //- Heat flux at surface field
114  volScalarField q_;
115 
116  //- Force density at surface field
117  volVectorField fD_;
118 
119  //- Number density field
120  volScalarField rhoN_;
121 
122  //- Mass density field
123  volScalarField rhoM_;
124 
125  //- Dsmc particle density field
126  volScalarField dsmcRhoN_;
127 
128  //- Linear kinetic energy density field
129  volScalarField linearKE_;
130 
131  //- Internal energy density field
132  volScalarField internalE_;
133 
134  // Internal degree of freedom density field
135  volScalarField iDof_;
136 
137  //- Momentum density field
138  volVectorField momentum_;
139 
140  //- Parcel constant properties - one for each type
142 
143  //- Random number generator
144  Random rndGen_;
145 
146 
147  // boundary value fields
148 
149  //- Boundary temperature
150  volScalarField boundaryT_;
151 
152  //- Boundary velocity
153  volVectorField boundaryU_;
154 
155 
156  // References to the cloud sub-models
157 
158  //- Binary collision model
160  binaryCollisionModel_;
161 
162  //- Wall interaction model
164  wallInteractionModel_;
165 
166  //- Inflow boundary model
168  inflowBoundaryModel_;
169 
170 
171  // Private Member Functions
172 
173  //- Build the constant properties for all of the species
174  void buildConstProps();
175 
176  //- Record which particles are in which cell
177  void buildCellOccupancy();
178 
179  //- Initialise the system
180  void initialise(const IOdictionary& dsmcInitialiseDict);
181 
182  //- Calculate collisions between molecules
183  void collisions();
184 
185  //- Reset the data accumulation field values to zero
186  void resetFields();
187 
188  //- Calculate the volume field data
189  void calculateFields();
190 
191 
192 public:
193 
194  // Constructors
195 
196  //- Construct given name and mesh, will read Parcels and fields from
197  // file
198  DSMCCloud
199  (
200  const word& cloudName,
201  const fvMesh& mesh,
202  bool readFields = true
203  );
204 
205  //- Construct given name, mesh and initialisation dictionary.
206  DSMCCloud
207  (
208  const word& cloudName,
209  const fvMesh& mesh,
210  const IOdictionary& dsmcInitialiseDict
211  );
212 
213  //- Disallow default bitwise copy construction
214  DSMCCloud(const DSMCCloud&) = delete;
215 
216 
217  //- Destructor
218  virtual ~DSMCCloud();
219 
220 
221  //- Type of parcel the cloud was instantiated for
222  typedef ParcelType parcelType;
223 
224 
225  // Member Functions
226 
227  // Access
228 
229  // References to the mesh and databases
230 
231  //- Return the cloud type
232  inline const word& cloudName() const;
233 
234  //- Return references to the mesh
235  inline const fvMesh& mesh() const;
236 
237 
238  // References to the dsmc specific data
239 
240  //- Return particle properties dictionary
241  inline const IOdictionary& particleProperties() const;
242 
243  //- Return the idList
244  inline const List<word>& typeIdList() const;
245 
246  //- Return the number of real particles represented by one
247  // parcel
248  inline scalar nParticle() const;
249 
250  //- Return the cell occupancy addressing
251  inline const List<DynamicList<ParcelType*>>&
252  cellOccupancy() const;
253 
254  //- Return the sigmaTcRMax field. non-const access to allow
255  // updating.
256  inline volScalarField& sigmaTcRMax();
257 
258  //- Return the collision selection remainder field. non-const
259  // access to allow updating.
261 
262  //- Return all of the constant properties
264  constProps() const;
265 
266  //- Return the constant properties of the given typeId
267  inline const typename ParcelType::constantProperties&
268  constProps(label typeId) const;
269 
270  //- Return references to the random object
271  inline Random& rndGen();
272 
273 
274  // References to the boundary fields for surface data collection
275 
276  //- Return non-const heat flux boundary field reference
277  inline volScalarField::Boundary& qBF();
278 
279  //- Return non-const force density at boundary field reference
280  inline volVectorField::Boundary& fDBF();
281 
282  //- Return non-const number density boundary field reference
283  inline volScalarField::Boundary& rhoNBF();
284 
285  //- Return non-const mass density boundary field reference
286  inline volScalarField::Boundary& rhoMBF();
287 
288  //- Return non-const linear kinetic energy density boundary
289  // field reference
290  inline volScalarField::Boundary& linearKEBF();
291 
292  //- Return non-const internal energy density boundary field
293  // reference
294  inline volScalarField::Boundary& internalEBF();
295 
296  //- Return non-const internal degree of freedom density boundary
297  // field reference
298  inline volScalarField::Boundary& iDofBF();
299 
300  //- Return non-const momentum density boundary field reference
301  inline volVectorField::Boundary& momentumBF();
302 
303 
304  // References to the macroscopic fields
305 
306  //- Return macroscopic temperature
307  inline const volScalarField& boundaryT() const;
308 
309  //- Return macroscopic velocity
310  inline const volVectorField& boundaryU() const;
311 
312  //- Return heat flux at surface field
313  inline const volScalarField& q() const;
314 
315  //- Return force density at surface field
316  inline const volVectorField& fD() const;
317 
318  //- Return the real particle number density field
319  inline const volScalarField& rhoN() const;
320 
321  //- Return the particle mass density field
322  inline const volScalarField& rhoM() const;
323 
324  //- Return the field of number of DSMC particles
325  inline const volScalarField& dsmcRhoN() const;
326 
327  //- Return the total linear kinetic energy (translational and
328  // thermal density field
329  inline const volScalarField& linearKE() const;
330 
331  //- Return the internal energy density field
332  inline const volScalarField& internalE() const;
333 
334  //- Return the average internal degrees of freedom field
335  inline const volScalarField& iDof() const;
336 
337  //- Return the momentum density field
338  inline const volVectorField& momentum() const;
339 
340 
341  // Kinetic theory helper functions
342 
343  //- Generate a random velocity sampled from the Maxwellian speed
344  // distribution
346  (
347  scalar temperature,
348  scalar mass
349  );
350 
351  //- Generate a random internal energy, sampled from the
352  // equilibrium distribution (Bird eqn 11.22 and 11.23 and
353  // adapting code from DSMC3.FOR)
355  (
356  scalar temperature,
357  direction internalDegreesOfFreedom
358  );
359 
360 
361  // From the Maxwellian distribution:
362  //- Average particle speed
363  inline scalar maxwellianAverageSpeed
364  (
365  scalar temperature,
366  scalar mass
367  ) const;
368 
370  (
371  scalarField temperature,
372  scalar mass
373  ) const;
374 
375  //- RMS particle speed
376  inline scalar maxwellianRMSSpeed
377  (
378  scalar temperature,
379  scalar mass
380  ) const;
381 
383  (
384  scalarField temperature,
385  scalar mass
386  ) const;
387 
388  //- Most probable speed
389  inline scalar maxwellianMostProbableSpeed
390  (
391  scalar temperature,
392  scalar mass
393  ) const;
394 
396  (
397  scalarField temperature,
398  scalar mass
399  ) const;
400 
401 
402  // Sub-models
403 
404  //- Return reference to binary elastic collision model
406  binaryCollision() const;
407 
408  //- Return non-const reference to binary elastic collision model
410  binaryCollision();
411 
412  //- Return reference to wall interaction model
414  wallInteraction() const;
415 
416  //- Return non-const reference to wall interaction model
418  wallInteraction();
419 
420  //- Return reference to wall interaction model
422  inflowBoundary() const;
423 
424  //- Return non-const reference to wall interaction model
426  inflowBoundary();
427 
428 
429  // Check
430 
431  //- Total mass in system
432  inline scalar massInSystem() const;
433 
434  //- Total linear momentum of the system
435  inline vector linearMomentumOfSystem() const;
436 
437  //- Total linear kinetic energy in the system
438  inline scalar linearKineticEnergyOfSystem() const;
439 
440  //- Total internal energy in the system
441  inline scalar internalEnergyOfSystem() const;
442 
443  //- Print cloud information
444  void info() const;
445 
446  //- Dump particle positions to .obj file
447  void dumpParticlePositions() const;
448 
449 
450 
451 
452  // Cloud evolution functions
453 
454  //- Add new parcel
455  void addNewParcel
456  (
457  const vector& position,
458  const label celli,
459  const vector& U,
460  const scalar Ei,
461  const label typeId
462  );
463 
464  //- Evolve the cloud (move, collide)
465  void evolve();
466 
467  //- Clear the Cloud
468  inline void clear();
469 
470 
471  // Mapping
472 
473  //- Remap the particles to the correct cells following mesh change
474  virtual void autoMap(const mapPolyMesh&);
475 
476 
477  // Member Operators
478 
479  //- Disallow default bitwise assignment
480  void operator=(const DSMCCloud&) = delete;
481 };
482 
483 
484 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
485 
486 } // End namespace Foam
487 
488 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
489 
490 #include "DSMCCloudI.H"
491 
492 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
493 
494 #ifdef NoRepository
495  #include "DSMCCloud.C"
496 #endif
497 
498 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
499 
500 #endif
501 
502 // ************************************************************************* //
TemplateName(blendedSchemeBase)
const volVectorField & momentum() const
Return the momentum density field.
Definition: DSMCCloudI.H:473
const volScalarField & boundaryT() const
Return macroscopic temperature.
Definition: DSMCCloudI.H:193
virtual void autoMap(const mapPolyMesh &)
Remap the particles to the correct cells following mesh change.
Definition: DSMCCloud.C:1087
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
const volScalarField & linearKE() const
Return the total linear kinetic energy (translational and.
Definition: DSMCCloudI.H:450
const volScalarField & dsmcRhoN() const
Return the field of number of DSMC particles.
Definition: DSMCCloudI.H:442
const List< word > & typeIdList() const
Return the idList.
Definition: DSMCCloudI.H:57
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
void info() const
Print cloud information.
Definition: DSMCCloud.C:970
vector linearMomentumOfSystem() const
Total linear momentum of the system.
Definition: DSMCCloudI.H:277
const word & cloudName() const
Return the cloud type.
Definition: DSMCCloudI.H:34
Templated inflow boundary model class.
Definition: DSMCCloud.H:61
const BinaryCollisionModel< DSMCCloud< ParcelType > > & binaryCollision() const
Return reference to binary elastic collision model.
Definition: DSMCCloudI.H:209
uint8_t direction
Definition: direction.H:45
volScalarField::Boundary & rhoMBF()
Return non-const mass density boundary field reference.
Definition: DSMCCloudI.H:153
Templated DSMC particle collision class.
Definition: DSMCCloud.H:55
scalar massInSystem() const
Total mass in system.
Definition: DSMCCloudI.H:256
volScalarField::Boundary & iDofBF()
Return non-const internal degree of freedom density boundary.
Definition: DSMCCloudI.H:177
const InflowBoundaryModel< DSMCCloud< ParcelType > > & inflowBoundary() const
Return reference to wall interaction model.
Definition: DSMCCloudI.H:241
scalar maxwellianRMSSpeed(scalar temperature, scalar mass) const
RMS particle speed.
Definition: DSMCCloudI.H:363
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
volScalarField::Boundary & qBF()
Return non-const heat flux boundary field reference.
Definition: DSMCCloudI.H:129
IOdictionary is derived from dictionary and IOobject to give the dictionary automatic IO functionalit...
Definition: IOdictionary.H:53
const volScalarField & internalE() const
Return the internal energy density field.
Definition: DSMCCloudI.H:458
const volVectorField & boundaryU() const
Return macroscopic velocity.
Definition: DSMCCloudI.H:201
A class for handling words, derived from string.
Definition: word.H:59
Random & rndGen()
Return references to the random object.
Definition: DSMCCloudI.H:121
const volScalarField & q() const
Return heat flux at surface field.
Definition: DSMCCloudI.H:412
volScalarField::Boundary & linearKEBF()
Return non-const linear kinetic energy density boundary.
Definition: DSMCCloudI.H:161
void operator=(const DSMCCloud &)=delete
Disallow default bitwise assignment.
volVectorField::Boundary & fDBF()
Return non-const force density at boundary field reference.
Definition: DSMCCloudI.H:137
scalar maxwellianMostProbableSpeed(scalar temperature, scalar mass) const
Most probable speed.
Definition: DSMCCloudI.H:388
volScalarField::Boundary & internalEBF()
Return non-const internal energy density boundary field.
Definition: DSMCCloudI.H:169
Random number generator.
Definition: Random.H:57
Base cloud calls templated on particle type.
Definition: Cloud.H:52
scalar equipartitionInternalEnergy(scalar temperature, direction internalDegreesOfFreedom)
Generate a random internal energy, sampled from the.
Definition: DSMCCloud.C:1027
const volScalarField & rhoM() const
Return the particle mass density field.
Definition: DSMCCloudI.H:434
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:221
void clear()
Clear the Cloud.
Definition: DSMCCloudI.H:480
const volScalarField & rhoN() const
Return the real particle number density field.
Definition: DSMCCloudI.H:427
scalar linearKineticEnergyOfSystem() const
Total linear kinetic energy in the system.
Definition: DSMCCloudI.H:299
scalar nParticle() const
Return the number of real particles represented by one.
Definition: DSMCCloudI.H:64
const volScalarField & iDof() const
Return the average internal degrees of freedom field.
Definition: DSMCCloudI.H:466
DSMCCloud(const word &cloudName, const fvMesh &mesh, bool readFields=true)
Construct given name and mesh, will read Parcels and fields from.
Definition: DSMCCloud.C:470
void addNewParcel(const vector &position, const label celli, const vector &U, const scalar Ei, const label typeId)
Add new parcel.
Definition: DSMCCloud.C:454
virtual ~DSMCCloud()
Destructor.
Definition: DSMCCloud.C:933
const fvMesh & mesh() const
Return references to the mesh.
Definition: DSMCCloudI.H:41
const volVectorField & fD() const
Return force density at surface field.
Definition: DSMCCloudI.H:419
U
Definition: pEqn.H:72
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
void evolve()
Evolve the cloud (move, collide)
Definition: DSMCCloud.C:940
scalar internalEnergyOfSystem() const
Total internal energy in the system.
Definition: DSMCCloudI.H:321
const List< DynamicList< ParcelType * > > & cellOccupancy() const
Return the cell occupancy addressing.
Definition: DSMCCloudI.H:72
volScalarField::Boundary & rhoNBF()
Return non-const number density boundary field reference.
Definition: DSMCCloudI.H:145
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
const WallInteractionModel< DSMCCloud< ParcelType > > & wallInteraction() const
Return reference to wall interaction model.
Definition: DSMCCloudI.H:225
const IOdictionary & particleProperties() const
Return particle properties dictionary.
Definition: DSMCCloudI.H:49
volScalarField & sigmaTcRMax()
Return the sigmaTcRMax field. non-const access to allow.
Definition: DSMCCloudI.H:79
Templated wall interaction model class.
Definition: DSMCCloud.H:58
void dumpParticlePositions() const
Dump particle positions to .obj file.
Definition: DSMCCloud.C:1063
scalarField & collisionSelectionRemainder()
Return the collision selection remainder field. non-const.
Definition: DSMCCloudI.H:87
vector equipartitionLinearVelocity(scalar temperature, scalar mass)
Generate a random velocity sampled from the Maxwellian speed.
Definition: DSMCCloud.C:1014
scalar maxwellianAverageSpeed(scalar temperature, scalar mass) const
Average particle speed.
Definition: DSMCCloudI.H:338
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:75
const List< typename ParcelType::constantProperties > & constProps() const
Return all of the constant properties.
Definition: DSMCCloudI.H:95
Namespace for OpenFOAM.
volVectorField::Boundary & momentumBF()
Return non-const momentum density boundary field reference.
Definition: DSMCCloudI.H:185