ReactingCloud.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-2021 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::ReactingCloud
26 
27 Description
28  Templated base class for reacting cloud
29 
30  - Adds to thermodynamic cloud
31  - Variable composition (single phase)
32  - Phase change
33 
34 SourceFiles
35  ReactingCloudI.H
36  ReactingCloud.C
37 
38 \*---------------------------------------------------------------------------*/
39 
40 #ifndef ReactingCloud_H
41 #define ReactingCloud_H
42 
43 #include "volFieldsFwd.H"
44 #include "fvMatricesFwd.H"
45 #include "dimensionedTypes.H"
46 #include "fvMesh.H"
47 #include "fluidThermo.H"
48 #include "Cloud.H"
49 
50 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
51 
52 namespace Foam
53 {
54 
55 // Forward declaration of classes
56 
57 template<class CloudType>
58 class PhaseChangeModel;
59 
60 
61 /*---------------------------------------------------------------------------*\
62  Class ReactingCloudName Declaration
63 \*---------------------------------------------------------------------------*/
64 
66 
67 
68 /*---------------------------------------------------------------------------*\
69  Class ReactingCloud Declaration
70 \*---------------------------------------------------------------------------*/
71 
72 template<class CloudType>
73 class ReactingCloud
74 :
75  public CloudType,
76  public ReactingCloudName
77 {
78 public:
79 
80  // Public Typedefs
81 
82  //- Type of cloud this cloud was instantiated for
83  typedef CloudType cloudType;
84 
85  //- Type of parcel the cloud was instantiated for
86  typedef typename CloudType::particleType parcelType;
87 
88  //- Convenience typedef for this cloud type
90 
91 
92 private:
93 
94  // Private Data
95 
96  //- Cloud copy pointer
97  autoPtr<ReactingCloud<CloudType>> cloudCopyPtr_;
98 
99 
100 
101 protected:
102 
103  // Protected data
104 
105  //- Parcel constant properties
106  typename parcelType::constantProperties constProps_;
107 
108 
109  // References to the cloud sub-models
110 
111  //- Reacting phase change model
114 
115 
116  // Sources
117 
118  //- Mass transfer fields - one per carrier phase specie
120 
121 
122  // Protected Member Functions
123 
124  // New parcel helper functions
125 
126  //- Check that size of a composition field is valid
128  (
129  const scalarField& YSupplied,
130  const scalarField& Y,
131  const word& YName
132  );
133 
134 
135  // Initialisation
136 
137  //- Set cloud sub-models
138  void setModels();
139 
140 
141  // Cloud evolution functions
142 
143  //- Reset state of cloud
145 
146 
147 public:
148 
149  // Constructors
150 
151  //- Construct given carrier fields and thermo
153  (
154  const word& cloudName,
155  const volScalarField& rho,
156  const volVectorField& U,
157  const dimensionedVector& g,
158  const fluidThermo& carrierThermo,
159  const bool readFields = true
160  );
161 
162  //- Copy constructor with new name
164 
165  //- Copy constructor with new name - creates bare cloud
167  (
168  const fvMesh& mesh,
169  const word& name,
171  );
172 
173  //- Disallow default bitwise copy construction
174  ReactingCloud(const ReactingCloud&) = delete;
175 
176  //- Construct and return clone based on (this) with new name
177  virtual autoPtr<Cloud<parcelType>> clone(const word& name)
178  {
180  (
181  new ReactingCloud(*this, name)
182  );
183  }
184 
185  //- Construct and return bare clone based on (this) with new name
186  virtual autoPtr<Cloud<parcelType>> cloneBare(const word& name) const
187  {
189  (
190  new ReactingCloud(this->mesh(), name, *this)
191  );
192  }
193 
194 
195  //- Destructor
196  virtual ~ReactingCloud();
197 
198 
199  // Member Functions
200 
201  // Access
202 
203  //- Return a reference to the cloud copy
204  inline const ReactingCloud& cloudCopy() const;
205 
206  //- Return the constant properties
207  inline const typename parcelType::constantProperties&
208  constProps() const;
209 
210  //- Return access to the constant properties
211  inline typename parcelType::constantProperties& constProps();
212 
213 
214  // Sub-models
215 
216  //- Return const access to reacting phase change model
218  phaseChange() const;
219 
220  //- Return reference to reacting phase change model
222  phaseChange();
223 
224 
225  // Sources
226 
227  //- Mass
228 
229  //- Return reference to mass source for field i
231  rhoTrans(const label i);
232 
233  //- Return const access to mass source fields
235  rhoTrans() const;
236 
237  //- Return reference to mass source fields
239  rhoTrans();
240 
241  //- Return mass source term for specie i - specie eqn
242  inline tmp<fvScalarMatrix> SYi
243  (
244  const label i,
245  const volScalarField& Yi
246  ) const;
247 
248  //- Return tmp total mass source for carrier phase
249  // - fully explicit
250  inline tmp<volScalarField::Internal> Srho() const;
251 
252  //- Return total mass source term [kg/m^3/s]
254  (
255  const volScalarField& rho
256  ) const;
257 
258 
259  // Cloud evolution functions
260 
261  //- Set parcel thermo properties
263  (
264  parcelType& parcel,
265  const scalar lagrangianDt
266  );
267 
268  //- Check parcel properties
270  (
271  parcelType& parcel,
272  const scalar lagrangianDt,
273  const bool fullyDescribed
274  );
275 
276  //- Store the current cloud state
277  void storeState();
278 
279  //- Reset the current cloud to the previously stored state
280  void restoreState();
281 
282  //- Reset the cloud source terms
283  void resetSourceTerms();
284 
285  //- Apply relaxation to (steady state) cloud sources
286  void relaxSources(const ReactingCloud<CloudType>& cloudOldTime);
287 
288  //- Apply scaling to (transient) cloud sources
289  void scaleSources();
290 
291  //- Evolve the cloud
292  void evolve();
293 
294 
295  // Mapping
296 
297  //- Remap the cells of particles corresponding to the
298  // mesh topology change with a default tracking data object
299  virtual void autoMap(const mapPolyMesh&);
300 
301 
302  // I-O
303 
304  //- Print cloud information
305  void info();
306 
307 
308  // Member Operators
309 
310  //- Disallow default bitwise assignment
311  void operator=(const ReactingCloud&) = delete;
312 };
313 
314 
315 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
316 
317 } // End namespace Foam
318 
319 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
320 
321 #include "ReactingCloudI.H"
322 
323 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
324 
325 #ifdef NoRepository
326  #include "ReactingCloud.C"
327 #endif
328 
329 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
330 
331 #endif
332 
333 // ************************************************************************* //
ParcelType particleType
Definition: Cloud.H:105
TemplateName(blendedSchemeBase)
virtual ~ReactingCloud()
Destructor.
void resetSourceTerms()
Reset the cloud source terms.
autoPtr< IOobject > clone() const
Clone.
Definition: IOobject.H:276
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 word & name() const
Return name.
Definition: IOobject.H:303
Templated phase change model class.
Definition: ReactingCloud.H:57
void storeState()
Store the current cloud state.
parcelType::constantProperties constProps_
Parcel constant properties.
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
const word & cloudName() const
Return the cloud type.
Definition: DSMCCloudI.H:34
const PtrList< volScalarField::Internal > & rhoTrans() const
Return const access to mass source fields.
void setParcelThermoProperties(parcelType &parcel, const scalar lagrangianDt)
Set parcel thermo properties.
void restoreState()
Reset the current cloud to the previously stored state.
tmp< fvScalarMatrix > SYi(const label i, const volScalarField &Yi) const
Return mass source term for specie i - specie eqn.
void setModels()
Set cloud sub-models.
Definition: ReactingCloud.C:33
tmp< volScalarField::Internal > Srho() const
Return tmp total mass source for carrier phase.
const dimensionedScalar c
Speed of light in a vacuum.
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:158
void evolve()
Evolve the cloud.
const ReactingCloud & cloudCopy() const
Return a reference to the cloud copy.
virtual autoPtr< Cloud< parcelType > > cloneBare(const word &name) const
Construct and return bare clone based on (this) with new name.
A class for handling words, derived from string.
Definition: word.H:59
Base-class for fluid thermodynamic properties.
Definition: fluidThermo.H:53
autoPtr< PhaseChangeModel< ReactingCloud< CloudType > > > phaseChangeModel_
Reacting phase change model.
const parcelType::constantProperties & constProps() const
Return the constant properties.
void checkParcelProperties(parcelType &parcel, const scalar lagrangianDt, const bool fullyDescribed)
Check parcel properties.
CloudType cloudType
Type of cloud this cloud was instantiated for.
Definition: ReactingCloud.H:82
virtual void autoMap(const mapPolyMesh &)
Remap the cells of particles corresponding to the.
void scaleSources()
Apply scaling to (transient) cloud sources.
const fvMesh & mesh() const
Return references to the mesh.
Definition: DSMCCloudI.H:41
Templated base class for reacting cloud.
Definition: ReactingCloud.H:72
Forward declarations of fvMatrix specialisations.
U
Definition: pEqn.H:72
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: List.H:70
PtrList< volScalarField > & Y
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:78
void info()
Print cloud information.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
void cloudReset(ReactingCloud< CloudType > &c)
Reset state of cloud.
Definition: ReactingCloud.C:67
PtrList< volScalarField::Internal > rhoTrans_
Mass transfer fields - one per carrier phase specie.
void checkSuppliedComposition(const scalarField &YSupplied, const scalarField &Y, const word &YName)
Check that size of a composition field is valid.
Definition: ReactingCloud.C:48
void operator=(const ReactingCloud &)=delete
Disallow default bitwise assignment.
An auto-pointer similar to the STL auto_ptr but with automatic casting to a reference to the type and...
Definition: PtrList.H:52
ReactingCloud(const word &cloudName, const volScalarField &rho, const volVectorField &U, const dimensionedVector &g, const fluidThermo &carrierThermo, const bool readFields=true)
Construct given carrier fields and thermo.
Definition: ReactingCloud.C:79
A class for managing temporary objects.
Definition: PtrList.H:53
const dimensionedVector & g
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:75
const PhaseChangeModel< ReactingCloud< CloudType > > & phaseChange() const
Return const access to reacting phase change model.
void relaxSources(const ReactingCloud< CloudType > &cloudOldTime)
Apply relaxation to (steady state) cloud sources.
CloudType::particleType parcelType
Type of parcel the cloud was instantiated for.
Definition: ReactingCloud.H:85
ReactingCloud< CloudType > reactingCloudType
Convenience typedef for this cloud type.
Definition: ReactingCloud.H:88
Namespace for OpenFOAM.