fractal.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) 2019-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::diameterModels::shapeModels::fractal
26 
27 Description
28  Class for modeling fractal shapes (e.g. of aerosol agglomerates) based on a
29  constant fractal dimension and an average but field-dependent surface area
30  to volume ratio kappa.
31 
32  The effect of sintering (coalescence of primary particles) on the surface
33  area is taken into account by a separate source term.
34 
35  Kappa is transported between the size groups as a secondary property by
36  means of coalescence (coagulation), breakup as well as drift.
37 
38  By assuming a monodisperse size distribution of the primary particles in the
39  aggregate, the collisional diameter of a size group can then be computed by
40 
41  \f[
42  d_{Coll_i} =
43  \frac{6}{\kappa_i}
44  \left(
45  \frac{v_i \kappa_i^3}{36 \pi \alpha_c}
46  \right)^{1/D_f}\;.
47  \f]
48 
49 Usage
50  \table
51  Property | Description | Required | Default value
52  kappa | Field and BC value | yes |
53  Df | Fractal dimension | yes |
54  alphaC | Constant | yes |
55  \endtable
56 
57 SourceFiles
58  fractal.C
59 
60 \*---------------------------------------------------------------------------*/
61 
62 #ifndef fractal_H
63 #define fractal_H
64 
65 #include "SecondaryPropertyModel.H"
66 #include "shapeModel.H"
67 
68 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
69 
70 namespace Foam
71 {
72 namespace diameterModels
73 {
74 namespace shapeModels
75 {
76 
77 class sinteringModel;
78 
79 /*---------------------------------------------------------------------------*\
80  Class fractal Declaration
81 \*---------------------------------------------------------------------------*/
82 
83 class fractal
84 :
85  public SecondaryPropertyModel<shapeModel>
86 {
87 public:
88 
89  // Public Data types
90 
91  //- Surface growth type enumeration
93  {
97  };
98 
99  //- Surface growth type names
101 
103 private:
104 
105  // Private Data
106 
107  //- Ratio of surface area to volume
108  volScalarField kappa_;
109 
110  //- Fractal dimension
112 
113  //- Scaling prefactor
116  //- Collisional diameter
117  volScalarField dColl_;
118 
119  //- Explicit source
120  volScalarField Su_;
121 
122  //- Sintering model
123  autoPtr<sinteringModel> sinteringModel_;
124 
125 
126  // Private Member Functions
127 
128  tmp<volScalarField> dColl() const;
129 
130 
131 public:
132 
133  //- Runtime type information
134  TypeName("fractal");
135 
136 
137  // Constructors
138 
139  //- Construct from dictionary and sizeGroup
140  fractal
141  (
142  const dictionary& dict,
143  const sizeGroup& group
144  );
145 
146  //- Disallow default bitwise copy construction
147  fractal(const fractal&) = delete;
148 
149 
150  //- Destructor
151  virtual ~fractal();
152 
153 
154  // Member Functions
155 
156  // Access
157 
158  //- Return reference to secondary property field
159  virtual const volScalarField& fld() const;
160 
161  //- Access to secondary property source
162  virtual volScalarField& src();
163 
164  //- Return representative surface area of the sizeGroup
165  virtual const tmp<volScalarField> a() const;
166 
167  //- Return representative diameter of the sizeGroup
168  virtual const tmp<volScalarField> d() const
169  {
170  return dColl_;
171  }
172 
173  // Edit
174 
175  //- Correct the collisional diameter
176  virtual void correct();
177 
178  //- Add drift contribution to secondary property source
179  virtual void addDrift
180  (
181  const volScalarField& Su,
182  const sizeGroup& fu,
183  const driftModel& model
184  );
185 
186 
187  // Member Operators
188 
189  //- Disallow default bitwise assignment
190  void operator=(const fractal&) = delete;
191 };
192 
193 
194 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
195 
196 } // End namespace shapeModels
197 } // End namespace diameterModels
198 } // End namespace Foam
199 
200 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
201 
202 #endif
203 
204 // ************************************************************************* //
surfaceGrowthTypes
Surface growth type enumeration.
Definition: fractal.H:111
dictionary dict
void operator=(const fractal &)=delete
Disallow default bitwise assignment.
Base class for drift models.
Definition: driftModel.H:50
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
fractal(const dictionary &dict, const sizeGroup &group)
Construct from dictionary and sizeGroup.
virtual void correct()
Correct the collisional diameter.
virtual const volScalarField & fld() const
Return reference to secondary property field.
word group() const
Return group (extension part of name)
Definition: IOobject.C:330
static const NamedEnum< surfaceGrowthTypes, 3 > sgTypeNames_
Surface growth type names.
Definition: fractal.H:119
virtual const tmp< volScalarField > d() const
Return representative diameter of the sizeGroup.
Definition: fractal.H:187
TypeName("fractal")
Runtime type information.
This class represents a single sizeGroup belonging to a velocityGroup. The main property of a sizeGro...
Definition: sizeGroup.H:100
virtual const tmp< volScalarField > a() const
Return representative surface area of the sizeGroup.
virtual volScalarField & src()
Access to secondary property source.
Class for modeling fractal shapes (e.g. of aerosol agglomerates) based on a constant fractal dimensio...
Definition: fractal.H:102
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 tmp< volScalarField::Internal > & Su
Definition: alphaSuSp.H:6
A class for managing temporary objects.
Definition: PtrList.H:53
virtual void addDrift(const volScalarField &Su, const sizeGroup &fu, const driftModel &model)
Add drift contribution to secondary property source.
Namespace for OpenFOAM.