ManualInjection.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2014 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 "ManualInjection.H"
27 #include "mathematicalConstants.H"
28 #include "PackedBoolList.H"
29 
30 using namespace Foam::constant::mathematical;
31 
32 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
34 template<class CloudType>
36 (
37  const dictionary& dict,
38  CloudType& owner,
39  const word& modelName
40 )
41 :
42  InjectionModel<CloudType>(dict, owner, modelName, typeName),
43  positionsFile_(this->coeffDict().lookup("positionsFile")),
44  positions_
45  (
46  IOobject
47  (
48  positionsFile_,
49  owner.db().time().constant(),
50  owner.mesh(),
53  )
54  ),
55  diameters_(positions_.size()),
56  injectorCells_(positions_.size(), -1),
57  injectorTetFaces_(positions_.size(), -1),
58  injectorTetPts_(positions_.size(), -1),
59  U0_(this->coeffDict().lookup("U0")),
60  sizeDistribution_
61  (
63  (
64  this->coeffDict().subDict("sizeDistribution"),
65  owner.rndGen()
66  )
67  ),
68  ignoreOutOfBounds_
69  (
70  this->coeffDict().lookupOrDefault("ignoreOutOfBounds", false)
71  )
72 {
73  updateMesh();
74 
75  // Construct parcel diameters
76  forAll(diameters_, i)
77  {
78  diameters_[i] = sizeDistribution_->sample();
79  }
80 
81  // Determine volume of particles to inject
82  this->volumeTotal_ = sum(pow3(diameters_))*pi/6.0;
83 }
84 
85 
86 template<class CloudType>
88 (
90 )
91 :
93  positionsFile_(im.positionsFile_),
94  positions_(im.positions_),
95  diameters_(im.diameters_),
96  injectorCells_(im.injectorCells_),
97  injectorTetFaces_(im.injectorTetFaces_),
98  injectorTetPts_(im.injectorTetPts_),
99  U0_(im.U0_),
100  sizeDistribution_(im.sizeDistribution_().clone().ptr()),
101  ignoreOutOfBounds_(im.ignoreOutOfBounds_)
102 {}
103 
104 
105 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
106 
107 template<class CloudType>
109 {}
110 
111 
112 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
113 
114 template<class CloudType>
116 {
117  label nRejected = 0;
118 
119  PackedBoolList keep(positions_.size(), true);
120 
121  forAll(positions_, pI)
122  {
123  if
124  (
125  !this->findCellAtPosition
126  (
127  injectorCells_[pI],
128  injectorTetFaces_[pI],
129  injectorTetPts_[pI],
130  positions_[pI],
131  !ignoreOutOfBounds_
132  )
133  )
134  {
135  keep[pI] = false;
136  nRejected++;
137  }
138  }
139 
140  if (nRejected > 0)
141  {
142  inplaceSubset(keep, positions_);
143  inplaceSubset(keep, diameters_);
144  inplaceSubset(keep, injectorCells_);
145  inplaceSubset(keep, injectorTetFaces_);
146  inplaceSubset(keep, injectorTetPts_);
147 
148  Info<< " " << nRejected
149  << " particles ignored, out of bounds" << endl;
150  }
151 }
152 
153 
154 template<class CloudType>
156 {
157  // Injection is instantaneous - but allow for a finite interval to
158  // avoid numerical issues when interval is zero
159  return ROOTVSMALL;
160 }
161 
162 
163 template<class CloudType>
165 (
166  const scalar time0,
167  const scalar time1
168 )
169 {
170  if ((0.0 >= time0) && (0.0 < time1))
171  {
172  return positions_.size();
173  }
174  else
175  {
176  return 0;
177  }
178 }
179 
180 
181 template<class CloudType>
183 (
184  const scalar time0,
185  const scalar time1
186 )
187 {
188  // All parcels introduced at SOI
189  if ((0.0 >= time0) && (0.0 < time1))
190  {
191  return this->volumeTotal_;
192  }
193  else
194  {
195  return 0.0;
196  }
197 }
198 
199 
200 template<class CloudType>
202 (
203  const label parcelI,
204  const label,
205  const scalar,
206  vector& position,
207  label& cellOwner,
208  label& tetFaceI,
209  label& tetPtI
210 )
211 {
212  position = positions_[parcelI];
213  cellOwner = injectorCells_[parcelI];
214  tetFaceI = injectorTetFaces_[parcelI];
215  tetPtI = injectorTetPts_[parcelI];
216 }
217 
218 
219 template<class CloudType>
221 (
222  const label parcelI,
223  const label,
224  const scalar,
225  typename CloudType::parcelType& parcel
226 )
227 {
228  // set particle velocity
229  parcel.U() = U0_;
230 
231  // set particle diameter
232  parcel.d() = diameters_[parcelI];
233 }
234 
235 
236 template<class CloudType>
238 {
239  return false;
240 }
241 
242 
243 template<class CloudType>
245 {
246  return true;
247 }
248 
249 
250 // ************************************************************************* //
dimensionedScalar pow3(const dimensionedScalar &ds)
virtual bool validInjection(const label parcelI)
Return flag to identify whether or not injection of parcelI is.
virtual label parcelsToInject(const scalar time0, const scalar time1)
Number of parcels to introduce relative to SOI.
virtual ~ManualInjection()
Destructor.
A bit-packed bool list.
A class for handling words, derived from string.
Definition: word.H:59
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
messageStream Info
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
virtual void setProperties(const label parcelI, const label nParcels, const scalar time, typename CloudType::parcelType &parcel)
Set the parcel properties.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:217
dictionary dict
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
virtual scalar volumeToInject(const scalar time0, const scalar time1)
Volume of parcels to introduce relative to SOI.
#define forAll(list, i)
Definition: UList.H:421
virtual void setPositionAndCell(const label parcelI, const label nParcels, const scalar time, vector &position, label &cellOwner, label &tetFaceI, label &tetPtI)
Set the injection position and owner cell, tetFace and tetPt.
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
Templated injection model class.
dimensioned< Type > sum(const DimensionedField< Type, GeoMesh > &df)
static autoPtr< distributionModel > New(const dictionary &dict, cachedRandom &rndGen)
Selector.
void inplaceSubset(const UList< T > &select, const T &value, ListType &)
Inplace extract elements of List when select is a certain value.
virtual void updateMesh()
Set injector locations when mesh is updated.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:68
mathematical constants.
scalar timeEnd() const
Return the end-of-injection time.
Manual injection.
ManualInjection(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.