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