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-2017 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  }
262 }
263 
264 
265 template<class CloudType>
267 {
268  return this->SOI_ + duration_;
269 }
270 
271 
272 template<class CloudType>
274 (
275  const scalar time0,
276  const scalar time1
277 )
278 {
279  if ((time0 >= 0.0) && (time0 < duration_))
280  {
281  return floor((time1 - time0)*parcelsPerSecond_);
282  }
283  else
284  {
285  return 0;
286  }
287 }
288 
289 
290 template<class CloudType>
292 (
293  const scalar time0,
294  const scalar time1
295 )
296 {
297  if ((time0 >= 0.0) && (time0 < duration_))
298  {
299  return flowRateProfile_.integrate(time0, time1);
300  }
301  else
302  {
303  return 0.0;
304  }
305 }
306 
307 
308 template<class CloudType>
310 (
311  const label,
312  const label,
313  const scalar,
314  vector& position,
315  label& cellOwner,
316  label& tetFacei,
317  label& tetPti
318 )
319 {
320  cachedRandom& rndGen = this->owner().rndGen();
321 
322  scalar beta = mathematical::twoPi*rndGen.sample01<scalar>();
323  normal_ = tanVec1_*cos(beta) + tanVec2_*sin(beta);
324 
325  switch (injectionMethod_)
326  {
327  case imPoint:
328  {
329  position = position_;
330  cellOwner = injectorCell_;
331  tetFacei = tetFacei_;
332  tetPti = tetPti_;
333 
334  break;
335  }
336  case imDisc:
337  {
338  scalar frac = rndGen.sample01<scalar>();
339  scalar dr = outerDiameter_ - innerDiameter_;
340  scalar r = 0.5*(innerDiameter_ + frac*dr);
341  position = position_ + r*normal_;
342 
343  this->findCellAtPosition
344  (
345  cellOwner,
346  tetFacei,
347  tetPti,
348  position,
349  false
350  );
351  break;
352  }
353  default:
354  {
356  << "Unknown injectionMethod type" << nl
357  << exit(FatalError);
358  }
359  }
360 }
361 
362 
363 template<class CloudType>
365 (
366  const label parcelI,
367  const label,
368  const scalar time,
369  typename CloudType::parcelType& parcel
370 )
371 {
372  cachedRandom& rndGen = this->owner().rndGen();
373 
374  // set particle velocity
375  const scalar deg2Rad = mathematical::pi/180.0;
376 
377  scalar t = time - this->SOI_;
378  scalar ti = thetaInner_.value(t);
379  scalar to = thetaOuter_.value(t);
380  scalar coneAngle = rndGen.sample01<scalar>()*(to - ti) + ti;
381 
382  coneAngle *= deg2Rad;
383  scalar alpha = sin(coneAngle);
384  scalar dcorr = cos(coneAngle);
385 
386  vector normal = alpha*normal_;
387  vector dirVec = dcorr*direction_;
388  dirVec += normal;
389  dirVec /= mag(dirVec);
390 
391  switch (flowType_)
392  {
393  case ftConstantVelocity:
394  {
395  parcel.U() = UMag_*dirVec;
396  break;
397  }
398  case ftPressureDrivenVelocity:
399  {
400  scalar pAmbient = this->owner().pAmbient();
401  scalar rho = parcel.rho();
402  scalar UMag = ::sqrt(2.0*(Pinj_.value(t) - pAmbient)/rho);
403  parcel.U() = UMag*dirVec;
404  break;
405  }
406  case ftFlowRateAndDischarge:
407  {
408  scalar Ao = 0.25*mathematical::pi*outerDiameter_*outerDiameter_;
409  scalar Ai = 0.25*mathematical::pi*innerDiameter_*innerDiameter_;
410  scalar massFlowRate =
411  this->massTotal()
412  *flowRateProfile_.value(t)
413  /this->volumeTotal();
414 
415  scalar Umag = massFlowRate/(parcel.rho()*Cd_.value(t)*(Ao - Ai));
416  parcel.U() = Umag*dirVec;
417  break;
418  }
419  default:
420  {
421  }
422  }
423 
424  // set particle diameter
425  parcel.d() = sizeDistribution_->sample();
426 }
427 
428 
429 template<class CloudType>
431 {
432  return false;
433 }
434 
435 
436 template<class CloudType>
438 {
439  return true;
440 }
441 
442 
443 // ************************************************************************* //
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
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)
Random number generator.
Definition: cachedRandom.H:63
virtual bool validInjection(const label parcelI)
Return flag to identify whether or not injection of parcelI is.
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
cachedRandom rndGen(label(0), -1)
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
static autoPtr< distributionModel > New(const dictionary &dict, cachedRandom &rndGen)
Selector.
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:262
scalar timeEnd() const
Return the end-of-injection time.
virtual void updateMesh()
Set injector locations when mesh is updated.
A normal distribution model.
virtual ~ConeNozzleInjection()
Destructor.
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