ConeNozzleInjection.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-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 \*---------------------------------------------------------------------------*/
25 
26 #include "ConeNozzleInjection.H"
27 #include "TimeFunction1.H"
28 #include "mathematicalConstants.H"
29 #include "distributionModel.H"
30 
31 using namespace Foam::constant;
32 
33 // * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * * //
34 
35 template<class CloudType>
37 {
38  word injectionMethodType = this->coeffDict().lookup("injectionMethod");
39  if (injectionMethodType == "disc")
40  {
41  injectionMethod_ = imDisc;
42  }
43  else if (injectionMethodType == "point")
44  {
45  injectionMethod_ = imPoint;
46 
47  // Set/cache the injector cell
48  this->findCellAtPosition
49  (
50  injectorCell_,
51  tetFacei_,
52  tetPti_,
53  position_,
54  false
55  );
56  }
57  else
58  {
60  << "injectionMethod must be either 'point' or 'disc'"
61  << exit(FatalError);
62  }
63 }
64 
65 
66 template<class CloudType>
68 {
69  word flowType = this->coeffDict().lookup("flowType");
70  if (flowType == "constantVelocity")
71  {
72  this->coeffDict().lookup("UMag") >> UMag_;
73  flowType_ = ftConstantVelocity;
74  }
75  else if (flowType == "pressureDrivenVelocity")
76  {
77  Pinj_.reset(this->coeffDict());
78  flowType_ = ftPressureDrivenVelocity;
79  }
80  else if (flowType == "flowRateAndDischarge")
81  {
82  Cd_.reset(this->coeffDict());
83  flowType_ = ftFlowRateAndDischarge;
84  }
85  else
86  {
88  << "flowType must be either 'constantVelocity', "
89  <<"'pressureDrivenVelocity' or 'flowRateAndDischarge'"
90  << exit(FatalError);
91  }
92 }
93 
94 
95 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
96 
97 template<class CloudType>
99 (
100  const dictionary& dict,
101  CloudType& owner,
102  const word& modelName
103 )
104 :
105  InjectionModel<CloudType>(dict, owner, modelName, typeName),
106  injectionMethod_(imPoint),
107  flowType_(ftConstantVelocity),
108  outerDiameter_(readScalar(this->coeffDict().lookup("outerDiameter"))),
109  innerDiameter_(readScalar(this->coeffDict().lookup("innerDiameter"))),
110  duration_(readScalar(this->coeffDict().lookup("duration"))),
111  position_(this->coeffDict().lookup("position")),
112  injectorCell_(-1),
113  tetFacei_(-1),
114  tetPti_(-1),
115  direction_(this->coeffDict().lookup("direction")),
116  parcelsPerSecond_
117  (
118  readScalar(this->coeffDict().lookup("parcelsPerSecond"))
119  ),
120  flowRateProfile_
121  (
123  (
124  owner.db().time(),
125  "flowRateProfile",
126  this->coeffDict()
127  )
128  ),
129  thetaInner_
130  (
132  (
133  owner.db().time(),
134  "thetaInner",
135  this->coeffDict()
136  )
137  ),
138  thetaOuter_
139  (
141  (
142  owner.db().time(),
143  "thetaOuter",
144  this->coeffDict()
145  )
146  ),
147  sizeDistribution_
148  (
150  (
151  this->coeffDict().subDict("sizeDistribution"),
152  owner.rndGen()
153  )
154  ),
155  tanVec1_(Zero),
156  tanVec2_(Zero),
157  normal_(Zero),
158 
159  UMag_(0.0),
160  Cd_(owner.db().time(), "Cd"),
161  Pinj_(owner.db().time(), "Pinj")
162 {
163  if (innerDiameter_ >= outerDiameter_)
164  {
166  << "innerNozzleDiameter >= outerNozzleDiameter" << nl
167  << exit(FatalError);
168  }
169 
170  duration_ = owner.db().time().userTimeToTime(duration_);
171 
172  setInjectionMethod();
173 
174  setFlowType();
175 
176  // Normalise direction vector
177  direction_ /= mag(direction_);
178 
179  // Determine direction vectors tangential to direction
180  tanVec1_ = normalised(perpendicular(direction_));
181  tanVec2_ = normalised(direction_ ^ tanVec1_);
182 
183  // Set total volume to inject
184  this->volumeTotal_ = flowRateProfile_.integrate(0.0, duration_);
185 
186  updateMesh();
187 }
188 
189 
190 template<class CloudType>
192 (
194 )
195 :
197  injectionMethod_(im.injectionMethod_),
198  flowType_(im.flowType_),
199  outerDiameter_(im.outerDiameter_),
200  innerDiameter_(im.innerDiameter_),
201  duration_(im.duration_),
202  position_(im.position_),
203  injectorCell_(im.injectorCell_),
204  tetFacei_(im.tetFacei_),
205  tetPti_(im.tetPti_),
206  direction_(im.direction_),
207  parcelsPerSecond_(im.parcelsPerSecond_),
208  flowRateProfile_(im.flowRateProfile_),
209  thetaInner_(im.thetaInner_),
210  thetaOuter_(im.thetaOuter_),
211  sizeDistribution_(im.sizeDistribution_().clone().ptr()),
212  tanVec1_(im.tanVec1_),
213  tanVec2_(im.tanVec2_),
214  normal_(im.normal_),
215  UMag_(im.UMag_),
216  Cd_(im.Cd_),
217  Pinj_(im.Pinj_)
218 {}
219 
220 
221 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
222 
223 template<class CloudType>
225 {}
226 
227 
228 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
229 
230 template<class CloudType>
232 {
233  // Set/cache the injector cells
234  switch (injectionMethod_)
235  {
236  case imPoint:
237  {
238  this->findCellAtPosition
239  (
240  injectorCell_,
241  tetFacei_,
242  tetPti_,
243  position_
244  );
245  }
246  default:
247  {}
248  }
249 }
250 
251 
252 template<class CloudType>
254 {
255  return this->SOI_ + duration_;
256 }
257 
258 
259 template<class CloudType>
261 (
262  const scalar time0,
263  const scalar time1
264 )
265 {
266  if ((time0 >= 0.0) && (time0 < duration_))
267  {
268  return floor((time1 - time0)*parcelsPerSecond_);
269  }
270  else
271  {
272  return 0;
273  }
274 }
275 
276 
277 template<class CloudType>
279 (
280  const scalar time0,
281  const scalar time1
282 )
283 {
284  if ((time0 >= 0.0) && (time0 < duration_))
285  {
286  return flowRateProfile_.integrate(time0, time1);
287  }
288  else
289  {
290  return 0.0;
291  }
292 }
293 
294 
295 template<class CloudType>
297 (
298  const label,
299  const label,
300  const scalar,
301  vector& position,
302  label& cellOwner,
303  label& tetFacei,
304  label& tetPti
305 )
306 {
307  Random& rndGen = this->owner().rndGen();
308 
309  scalar beta = mathematical::twoPi*rndGen.globalScalar01();
310  normal_ = tanVec1_*cos(beta) + tanVec2_*sin(beta);
311 
312  switch (injectionMethod_)
313  {
314  case imPoint:
315  {
316  position = position_;
317  cellOwner = injectorCell_;
318  tetFacei = tetFacei_;
319  tetPti = tetPti_;
320 
321  break;
322  }
323  case imDisc:
324  {
325  scalar frac = rndGen.globalScalar01();
326  scalar dr = outerDiameter_ - innerDiameter_;
327  scalar r = 0.5*(innerDiameter_ + frac*dr);
328  position = position_ + r*normal_;
329 
330  this->findCellAtPosition
331  (
332  cellOwner,
333  tetFacei,
334  tetPti,
335  position,
336  false
337  );
338  break;
339  }
340  default:
341  {
343  << "Unknown injectionMethod type" << nl
344  << exit(FatalError);
345  }
346  }
347 }
348 
349 
350 template<class CloudType>
352 (
353  const label parcelI,
354  const label,
355  const scalar time,
356  typename CloudType::parcelType& parcel
357 )
358 {
359  Random& rndGen = this->owner().rndGen();
360 
361  // set particle velocity
362  const scalar deg2Rad = mathematical::pi/180.0;
363 
364  scalar t = time - this->SOI_;
365  scalar ti = thetaInner_.value(t);
366  scalar to = thetaOuter_.value(t);
367  scalar coneAngle = rndGen.sample01<scalar>()*(to - ti) + ti;
368 
369  coneAngle *= deg2Rad;
370  scalar alpha = sin(coneAngle);
371  scalar dcorr = cos(coneAngle);
372 
373  vector normal = alpha*normal_;
374  vector dirVec = dcorr*direction_;
375  dirVec += normal;
376  dirVec /= mag(dirVec);
377 
378  switch (flowType_)
379  {
380  case ftConstantVelocity:
381  {
382  parcel.U() = UMag_*dirVec;
383  break;
384  }
385  case ftPressureDrivenVelocity:
386  {
387  scalar pAmbient = this->owner().pAmbient();
388  scalar rho = parcel.rho();
389  scalar UMag = ::sqrt(2.0*(Pinj_.value(t) - pAmbient)/rho);
390  parcel.U() = UMag*dirVec;
391  break;
392  }
393  case ftFlowRateAndDischarge:
394  {
395  scalar Ao = 0.25*mathematical::pi*outerDiameter_*outerDiameter_;
396  scalar Ai = 0.25*mathematical::pi*innerDiameter_*innerDiameter_;
397  scalar massFlowRate =
398  this->massTotal()
399  *flowRateProfile_.value(t)
400  /this->volumeTotal();
401 
402  scalar Umag = massFlowRate/(parcel.rho()*Cd_.value(t)*(Ao - Ai));
403  parcel.U() = Umag*dirVec;
404  break;
405  }
406  default:
407  {
408  }
409  }
410 
411  // set particle diameter
412  parcel.d() = sizeDistribution_->sample();
413 }
414 
415 
416 template<class CloudType>
418 {
419  return false;
420 }
421 
422 
423 template<class CloudType>
425 {
426  return true;
427 }
428 
429 
430 // ************************************************************************* //
Collection of constants.
virtual void setProperties(const label parcelI, const label nParcels, const scalar time, typename CloudType::parcelType &parcel)
Set the parcel properties.
dictionary dict
VectorSpace< Form, Cmpt, Ncmpts > normalised(const VectorSpace< Form, Cmpt, Ncmpts > &vs)
Definition: VectorSpaceI.H:413
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
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
error FatalError
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:319
Templated injection model class.
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.
dimensionedScalar sqrt(const dimensionedScalar &ds)
virtual bool validInjection(const label parcelI)
Return flag to identify whether or not injection of parcelI is.
Random rndGen(label(0))
virtual scalar volumeToInject(const scalar time0, const scalar time1)
Volume of parcels to introduce relative to SOI.
Type sample01()
Advance the state and return a sample of a given type from a.
Definition: RandomI.H:70
stressControl lookup("compactNormalStress") >> compactNormalStress
dimensionedScalar cos(const dimensionedScalar &ds)
virtual label parcelsToInject(const scalar time0, const scalar time1)
Number of parcels to introduce relative to SOI.
A class for handling words, derived from string.
Definition: word.H:59
scalar globalScalar01()
Advance the state and return a scalar sample from a uniform.
Definition: Random.C:94
static const zero Zero
Definition: zero.H:97
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:68
Random number generator.
Definition: Random.H:57
const scalar twoPi(2 *pi)
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:218
dimensionedScalar sin(const dimensionedScalar &ds)
static const char nl
Definition: Ostream.H:265
scalar timeEnd() const
Return the end-of-injection time.
virtual void updateMesh()
Set injector locations when mesh is updated.
Vector< Cmpt > perpendicular(const Vector< Cmpt > &v)
Definition: VectorI.H:166
virtual ~ConeNozzleInjection()
Destructor.
static autoPtr< distributionModel > New(const dictionary &dict, Random &rndGen)
Selector.
dimensionedScalar beta("beta", dimless/dimTemperature, laminarTransport)
dimensioned< scalar > mag(const dimensioned< Type > &)
ConeNozzleInjection(const dictionary &dict, CloudType &owner, const word &modelName)
Construct from dictionary.
const dimensionedScalar alpha
Fine-structure constant: default SI units: [].
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:69