ConeNozzleInjection.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-2016 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  cachedRandom& rndGen = this->owner().rndGen();
177 
178  // Normalise direction vector
179  direction_ /= mag(direction_);
180 
181  // Determine direction vectors tangential to direction
182  vector tangent = Zero;
183  scalar magTangent = 0.0;
184 
185  while(magTangent < SMALL)
186  {
187  vector v = rndGen.sample01<vector>();
188 
189  tangent = v - (v & direction_)*direction_;
190  magTangent = mag(tangent);
191  }
192 
193  tanVec1_ = tangent/magTangent;
194  tanVec2_ = direction_^tanVec1_;
195 
196  // Set total volume to inject
197  this->volumeTotal_ = flowRateProfile_.integrate(0.0, duration_);
198 
199  updateMesh();
200 }
201 
202 
203 template<class CloudType>
205 (
207 )
208 :
210  injectionMethod_(im.injectionMethod_),
211  flowType_(im.flowType_),
212  outerDiameter_(im.outerDiameter_),
213  innerDiameter_(im.innerDiameter_),
214  duration_(im.duration_),
215  position_(im.position_),
216  injectorCell_(im.injectorCell_),
217  tetFacei_(im.tetFacei_),
218  tetPtI_(im.tetPtI_),
219  direction_(im.direction_),
220  parcelsPerSecond_(im.parcelsPerSecond_),
221  flowRateProfile_(im.flowRateProfile_),
222  thetaInner_(im.thetaInner_),
223  thetaOuter_(im.thetaOuter_),
224  sizeDistribution_(im.sizeDistribution_().clone().ptr()),
225  tanVec1_(im.tanVec1_),
226  tanVec2_(im.tanVec2_),
227  normal_(im.normal_),
228  UMag_(im.UMag_),
229  Cd_(im.Cd_),
230  Pinj_(im.Pinj_)
231 {}
232 
233 
234 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
235 
236 template<class CloudType>
238 {}
239 
240 
241 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
242 
243 template<class CloudType>
245 {
246  // Set/cache the injector cells
247  switch (injectionMethod_)
248  {
249  case imPoint:
250  {
251  this->findCellAtPosition
252  (
253  injectorCell_,
254  tetFacei_,
255  tetPtI_,
256  position_
257  );
258  }
259  default:
260  {
261  // do nothing
262  }
263  }
264 }
265 
266 
267 template<class CloudType>
269 {
270  return this->SOI_ + duration_;
271 }
272 
273 
274 template<class CloudType>
276 (
277  const scalar time0,
278  const scalar time1
279 )
280 {
281  if ((time0 >= 0.0) && (time0 < duration_))
282  {
283  return floor((time1 - time0)*parcelsPerSecond_);
284  }
285  else
286  {
287  return 0;
288  }
289 }
290 
291 
292 template<class CloudType>
294 (
295  const scalar time0,
296  const scalar time1
297 )
298 {
299  if ((time0 >= 0.0) && (time0 < duration_))
300  {
301  return flowRateProfile_.integrate(time0, time1);
302  }
303  else
304  {
305  return 0.0;
306  }
307 }
308 
309 
310 template<class CloudType>
312 (
313  const label,
314  const label,
315  const scalar,
316  vector& position,
317  label& cellOwner,
318  label& tetFacei,
319  label& tetPtI
320 )
321 {
322  cachedRandom& rndGen = this->owner().rndGen();
323 
324  scalar beta = mathematical::twoPi*rndGen.sample01<scalar>();
325  normal_ = tanVec1_*cos(beta) + tanVec2_*sin(beta);
326 
327  switch (injectionMethod_)
328  {
329  case imPoint:
330  {
331  position = position_;
332  cellOwner = injectorCell_;
333  tetFacei = tetFacei_;
334  tetPtI = tetPtI_;
335 
336  break;
337  }
338  case imDisc:
339  {
340  scalar frac = rndGen.sample01<scalar>();
341  scalar dr = outerDiameter_ - innerDiameter_;
342  scalar r = 0.5*(innerDiameter_ + frac*dr);
343  position = position_ + r*normal_;
344 
345  this->findCellAtPosition
346  (
347  cellOwner,
348  tetFacei,
349  tetPtI,
350  position,
351  false
352  );
353  break;
354  }
355  default:
356  {
358  << "Unknown injectionMethod type" << nl
359  << exit(FatalError);
360  }
361  }
362 }
363 
364 
365 template<class CloudType>
367 (
368  const label parcelI,
369  const label,
370  const scalar time,
371  typename CloudType::parcelType& parcel
372 )
373 {
374  cachedRandom& rndGen = this->owner().rndGen();
375 
376  // set particle velocity
377  const scalar deg2Rad = mathematical::pi/180.0;
378 
379  scalar t = time - this->SOI_;
380  scalar ti = thetaInner_.value(t);
381  scalar to = thetaOuter_.value(t);
382  scalar coneAngle = rndGen.sample01<scalar>()*(to - ti) + ti;
383 
384  coneAngle *= deg2Rad;
385  scalar alpha = sin(coneAngle);
386  scalar dcorr = cos(coneAngle);
387 
388  vector normal = alpha*normal_;
389  vector dirVec = dcorr*direction_;
390  dirVec += normal;
391  dirVec /= mag(dirVec);
392 
393  switch (flowType_)
394  {
395  case ftConstantVelocity:
396  {
397  parcel.U() = UMag_*dirVec;
398  break;
399  }
400  case ftPressureDrivenVelocity:
401  {
402  scalar pAmbient = this->owner().pAmbient();
403  scalar rho = parcel.rho();
404  scalar UMag = ::sqrt(2.0*(Pinj_.value(t) - pAmbient)/rho);
405  parcel.U() = UMag*dirVec;
406  break;
407  }
408  case ftFlowRateAndDischarge:
409  {
410  scalar Ao = 0.25*mathematical::pi*outerDiameter_*outerDiameter_;
411  scalar Ai = 0.25*mathematical::pi*innerDiameter_*innerDiameter_;
412  scalar massFlowRate =
413  this->massTotal()
414  *flowRateProfile_.value(t)
415  /this->volumeTotal();
416 
417  scalar Umag = massFlowRate/(parcel.rho()*Cd_.value(t)*(Ao - Ai));
418  parcel.U() = Umag*dirVec;
419  break;
420  }
421  default:
422  {
423  }
424  }
425 
426  // set particle diameter
427  parcel.d() = sizeDistribution_->sample();
428 }
429 
430 
431 template<class CloudType>
433 {
434  return false;
435 }
436 
437 
438 template<class CloudType>
440 {
441  return true;
442 }
443 
444 
445 // ************************************************************************* //
Collection of constants.
virtual void setProperties(const label parcelI, const label nParcels, const scalar time, typename CloudType::parcelType &parcel)
Set the parcel properties.
cachedRandom rndGen(label(0),-1)
dictionary dict
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.
dimensionedScalar sqrt(const dimensionedScalar &ds)
Random number generator.
Definition: cachedRandom.H:63
virtual bool validInjection(const label parcelI)
Return flag to identify whether or not injection of parcelI is.
static autoPtr< distributionModel > New(const dictionary &dict, cachedRandom &rndGen)
Selector.
Type sample01()
Return a sample whose components lie in the range 0-1.
virtual scalar volumeToInject(const scalar time0, const scalar time1)
Volume of parcels to introduce relative to SOI.
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
static const zero Zero
Definition: zero.H:91
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if succesful.
Definition: doubleScalar.H:63
const scalar twoPi(2 *pi)
virtual bool fullyDescribed() const
Flag to identify whether model fully describes the parcel.
ParcelType parcelType
Type of parcel the cloud was instantiated for.
Definition: DSMCCloud.H:217
dimensionedScalar sin(const dimensionedScalar &ds)
static const char nl
Definition: Ostream.H:262
virtual void updateMesh()
Set injector locations when mesh is updated.
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.
A normal distribution model.
virtual ~ConeNozzleInjection()
Destructor.
scalar timeEnd() const
Return the end-of-injection time.
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: [].
Templated base class for dsmc cloud.
Definition: DSMCCloud.H:68