ParticleForceList.C
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-2024 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 "ParticleForceList.H"
27 #include "entry.H"
28 
29 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
30 
31 template<class CloudType>
33 (
34  CloudType& owner,
35  const fvMesh& mesh
36 )
37 :
39  owner_(owner),
40  mesh_(mesh),
41  dict_(dictionary::null),
42  calcCoupled_(true),
43  calcNonCoupled_(true)
44 {}
45 
46 
47 template<class CloudType>
49 (
50  CloudType& owner,
51  const fvMesh& mesh,
52  const dictionary& dict,
53  const bool readFields
54 )
55 :
57  owner_(owner),
58  mesh_(mesh),
59  dict_(dict),
60  calcCoupled_(true),
61  calcNonCoupled_(true)
62 {
63  if (readFields)
64  {
65  wordList modelNames(dict.toc());
66 
67  Info<< "Constructing particle forces" << endl;
68 
69  if (modelNames.size() > 0)
70  {
71  this->setSize(modelNames.size());
72 
73  label i = 0;
75  {
76  const word& model = iter().keyword();
77  if (iter().isDict())
78  {
79  this->set
80  (
81  i++,
83  (
84  owner,
85  mesh,
86  iter().dict(),
87  model
88  )
89  );
90  }
91  else
92  {
93  this->set
94  (
95  i++,
97  (
98  owner,
99  mesh,
101  model
102  )
103  );
104  }
105  }
106  }
107  else
108  {
109  Info<< " none" << endl;
110  }
111  }
112 }
113 
114 
115 template<class CloudType>
117 (
118  const ParticleForceList& pf
119 )
120 :
122  owner_(pf.owner_),
123  mesh_(pf.mesh_),
124  dict_(pf.dict_),
125  calcCoupled_(pf.calcCoupled_),
126  calcNonCoupled_(pf.calcNonCoupled_)
127 {}
128 
129 
130 // * * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * //
131 
132 template<class CloudType>
134 {}
135 
136 
137 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
138 
139 template<class CloudType>
141 {
142  forAll(*this, i)
143  {
144  this->operator[](i).cacheFields(store);
145  }
146 }
147 
148 
149 template<class CloudType>
151 (
152  const typename CloudType::parcelType& p,
153  const typename CloudType::parcelType::trackingData& td,
154  const scalar dt,
155  const scalar mass,
156  const scalar Re,
157  const scalar muc
158 ) const
159 {
160  forceSuSp value(Zero, 0.0);
161 
162  if (calcCoupled_)
163  {
164  forAll(*this, i)
165  {
166  value += this->operator[](i).calcCoupled(p, td, dt, mass, Re, muc);
167  }
168  }
169 
170  return value;
171 }
172 
173 
174 template<class CloudType>
176 (
177  const typename CloudType::parcelType& p,
178  const typename CloudType::parcelType::trackingData& td,
179  const scalar dt,
180  const scalar mass,
181  const scalar Re,
182  const scalar muc
183 ) const
184 {
185  forceSuSp value(Zero, 0.0);
186 
187  if (calcNonCoupled_)
188  {
189  forAll(*this, i)
190  {
191  value +=
192  this->operator[](i).calcNonCoupled(p, td, dt, mass, Re, muc);
193  }
194  }
195 
196  return value;
197 }
198 
199 
200 template<class CloudType>
202 (
203  const typename CloudType::parcelType& p,
204  const typename CloudType::parcelType::trackingData& td,
205  const scalar mass
206 ) const
207 {
208  scalar massEff = mass;
209  forAll(*this, i)
210  {
211  massEff += this->operator[](i).massAdd(p, td, mass);
212  }
213 
214  return massEff;
215 }
216 
217 
218 // ************************************************************************* //
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define forAllConstIter(Container, container, iter)
Iterate across all elements in the container object of type.
Definition: UList.H:477
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:80
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:225
Template class for intrusive linked lists.
Definition: ILList.H:67
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
List of particle forces.
ParticleForceList(CloudType &owner, const fvMesh &mesh)
Null constructor.
virtual forceSuSp calcNonCoupled(const typename CloudType::parcelType &p, const typename CloudType::parcelType::trackingData &td, const scalar dt, const scalar mass, const scalar Re, const scalar muc) const
Calculate the non-coupled force.
virtual scalar massEff(const typename CloudType::parcelType &p, const typename CloudType::parcelType::trackingData &td, const scalar mass) const
Return the effective mass.
virtual void cacheFields(const bool store)
Cache fields.
const dictionary & dict() const
Return the forces dictionary.
const CloudType & owner() const
Return const access to the cloud owner.
virtual ~ParticleForceList()
Destructor.
const fvMesh & mesh() const
Return the mesh database.
virtual forceSuSp calcCoupled(const typename CloudType::parcelType &p, const typename CloudType::parcelType::trackingData &td, const scalar dt, const scalar mass, const scalar Re, const scalar muc) const
Calculate the coupled force.
Abstract base class for particle forces.
Definition: ParticleForce.H:54
A templated 1D list of pointers to objects of type <T>, where the size of the array is known and used...
Definition: PtrList.H:75
bool set(const label) const
Is element set.
Definition: PtrListI.H:62
void setSize(const label)
Reset size of PtrList. If extending the PtrList, new entries are.
Definition: PtrList.C:131
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:162
wordList toc() const
Return the table of contents.
Definition: dictionary.C:976
static const dictionary null
Null dictionary.
Definition: dictionary.H:258
Helper container for force Su and Sp terms.
Definition: forceSuSp.H:64
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:99
A class for handling words, derived from string.
Definition: word.H:62
static const zero Zero
Definition: zero.H:97
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
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:257
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
messageStream Info
scalarField Re(const UList< complex > &cf)
Definition: complexFields.C:97
dictionary dict
volScalarField & p