FreeFOAM The Cross-Platform CFD Toolkit
multiHoleInjector.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) 1991-2010 OpenCFD Ltd.
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 "multiHoleInjector.H"
28 #include <OpenFOAM/Random.H>
30 
31 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
32 namespace Foam
33 {
34 
35 defineTypeNameAndDebug(multiHoleInjector, 0);
36 
38 (
39  injectorType,
40  multiHoleInjector,
41  dictionary
42 );
43 }
44 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
45 
46 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
47 
48 // Construct from components
49 Foam::multiHoleInjector::multiHoleInjector
50 (
51  const Foam::Time& t,
52  const Foam::dictionary& dict
53 )
54 :
55  injectorType(t, dict),
56  propsDict_(dict.subDict(typeName + "Props")),
57  centerPosition_(propsDict_.lookup("position")),
58  xyAngle_(readScalar(propsDict_.lookup("xyAngle"))),
59  zAngle_(readScalar(propsDict_.lookup("zAngle"))),
60  nHoles_(readLabel(propsDict_.lookup("nHoles"))),
61  umbrellaAngle_(readScalar(propsDict_.lookup("umbrellaAngle"))),
62  nozzleTipDiameter_(readScalar(propsDict_.lookup("nozzleTipDiameter"))),
63  angleSpacing_(propsDict_.lookup("angleSpacing")),
64  d_(readScalar(propsDict_.lookup("diameter"))),
65  Cd_(readScalar(propsDict_.lookup("Cd"))),
66  mass_(readScalar(propsDict_.lookup("mass"))),
67  nParcels_(readLabel(propsDict_.lookup("nParcels"))),
68  X_(propsDict_.lookup("X")),
69  massFlowRateProfile_(propsDict_.lookup("massFlowRateProfile")),
70  velocityProfile_(massFlowRateProfile_),
71  injectionPressureProfile_(massFlowRateProfile_),
72  CdProfile_(massFlowRateProfile_),
73  TProfile_(propsDict_.lookup("temperatureProfile")),
74  averageParcelMass_(nHoles_*mass_/nParcels_),
75  direction_(nHoles_),
76  position_(nHoles_),
77  pressureIndependentVelocity_(true),
78  tangentialInjectionVector1_(nHoles_),
79  tangentialInjectionVector2_(nHoles_)
80 {
81 
82 
83  // check if time entries for soi and eoi match
84  if (mag(massFlowRateProfile_[0][0]-TProfile_[0][0]) > SMALL)
85  {
86  FatalError << "multiHoleInjector::multiHoleInjector(const time& t, const dictionary dict) " << endl
87  << " start-times do not match for TemperatureProfile and massFlowRateProfile."
88  << abort(FatalError);
89  }
90 
91  if (mag(massFlowRateProfile_[massFlowRateProfile_.size()-1][0]-TProfile_[TProfile_.size()-1][0]) > SMALL)
92  {
93  FatalError << "multiHoleInjector::multiHoleInjector(const time& t, const dictionary dict) " << endl
94  << " end-times do not match for TemperatureProfile and massFlowRateProfile."
95  << abort(FatalError);
96  }
97 
98  // convert CA to real time
99  forAll(massFlowRateProfile_, i)
100  {
101  massFlowRateProfile_[i][0] = t.userTimeToTime(massFlowRateProfile_[i][0]);
102  velocityProfile_[i][0] = massFlowRateProfile_[i][0];
103  injectionPressureProfile_[i][0] = massFlowRateProfile_[i][0];
104  }
105 
106  forAll(TProfile_, i)
107  {
108  TProfile_[i][0] = t.userTimeToTime(TProfile_[i][0]);
109  }
110 
111  scalar integratedMFR = integrateTable(massFlowRateProfile_);
112 
113  forAll(massFlowRateProfile_, i)
114  {
115  // correct the massFlowRateProfile to match the injected mass
116  massFlowRateProfile_[i][1] *= mass_/integratedMFR;
117 
118  CdProfile_[i][0] = massFlowRateProfile_[i][0];
119  CdProfile_[i][1] = Cd_;
120  }
121 
122  setTangentialVectors();
123 
124  // check molar fractions
125  scalar Xsum = 0.0;
126  forAll(X_, i)
127  {
128  Xsum += X_[i];
129  }
130 
131  if (mag(Xsum - 1.0) > SMALL)
132  {
133  Info << "Warning!!!\n multiHoleInjector::multiHoleInjector(const time& t, Istream& is)"
134  << "X does not add up to 1.0, correcting molar fractions."
135  << endl;
136  forAll(X_, i)
137  {
138  X_[i] /= Xsum;
139  }
140  }
141 
142 }
143 
144 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
145 
147 {}
148 
149 
150 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
151 
152 void Foam::multiHoleInjector::setTangentialVectors()
153 {
154  scalar pi180 = mathematicalConstant::pi/180.0;
155  scalar alpha = xyAngle_*pi180;
156  scalar phi = zAngle_*pi180;
157 
158  vector xp(cos(alpha), sin(alpha), 0.0);
159  vector zp(cos(alpha)*sin(phi), sin(alpha)*sin(phi), cos(phi));
160  if (mag(zp-xp) < 1.0e-15)
161  {
162  xp = vector(0.0, 0.0, -1.0);
163  xp -= (xp & zp)*zp;
164  xp /= mag(xp);
165  }
166  vector yp = zp^xp;
167 
168 // Info << "xp = " << xp << endl;
169 // Info << "yp = " << yp << endl;
170 // Info << "zp = " << zp << endl;
171 
172  scalar angle = 0.0;
173  scalar u = umbrellaAngle_*pi180/2.0;
174  for(label i=0; i<nHoles_; i++)
175  {
176  angle += angleSpacing_[i];
177  scalar v = angle*pi180;
178  direction_[i] = cos(v)*sin(u)*xp + sin(v)*sin(u)*yp + cos(u)*zp;
179  vector dp = direction_[i] - (direction_[i] & zp)*direction_[i];
180  if (mag(dp) > SMALL)
181  {
182  dp /= mag(dp);
183  }
184  position_[i] = centerPosition_ + 0.5*nozzleTipDiameter_*dp;
185 // Info << "i = " << i << ", dir = " << direction_[i] << ", pos = " << position_[i] << endl;
186  }
187 
188  Random rndGen(label(0));
189 
190  for(label i=0; i<nHoles_; i++)
191  {
192  vector tangent(vector::zero);
193  scalar magV = 0;
194  while (magV < SMALL)
195  {
196  vector testThis = rndGen.vector01();
197 
198  tangent = testThis - (testThis & direction_[i])*direction_[i];
199  magV = mag(tangent);
200  }
201 
202  tangentialInjectionVector1_[i] = tangent/magV;
203  tangentialInjectionVector2_[i] = direction_[i] ^ tangentialInjectionVector1_[i];
204 
205  }
206 }
207 
208 
210 (
211  const scalar time0,
212  const scalar time1
213 ) const
214 {
215 
216  scalar mInj = mass_*(fractionOfInjection(time1)-fractionOfInjection(time0));
217  label nParcels = label(mInj/averageParcelMass_ + 0.49);
218 
219  return nParcels;
220 }
221 
223 {
224  return position_[n];
225 }
226 
228 (
229  const label n,
230  const scalar time,
231  const bool twoD,
232  const scalar angleOfWedge,
233  const vector& axisOfSymmetry,
234  const vector& axisOfWedge,
235  const vector& axisOfWedgeNormal,
236  Random& rndGen
237 ) const
238 {
239  if (twoD)
240  {
241  scalar is = position_[n] & axisOfSymmetry;
242  scalar magInj = mag(position_[n] - is*axisOfSymmetry);
243 
244  vector halfWedge =
245  axisOfWedge*cos(0.5*angleOfWedge)
246  + axisOfWedgeNormal*sin(0.5*angleOfWedge);
247  halfWedge /= mag(halfWedge);
248 
249  return (is*axisOfSymmetry + magInj*halfWedge);
250  }
251  else
252  {
253  // otherwise, disc injection
254  scalar iRadius = d_*rndGen.scalar01();
255  scalar iAngle = 2.0*mathematicalConstant::pi*rndGen.scalar01();
256 
257  return
258  (
259  position_[n]
260  + iRadius
261  * (
262  tangentialInjectionVector1_[n]*cos(iAngle)
263  + tangentialInjectionVector2_[n]*sin(iAngle)
264  )
265  );
266 
267  }
268  return position_[0];
269 }
270 
272 {
273  return nHoles_;
274 }
275 
276 Foam::scalar Foam::multiHoleInjector::d() const
277 {
278  return d_;
279 }
280 
282 (
283  const label i,
284  const scalar time
285 ) const
286 {
287  return direction_[i];
288 }
289 
291 (
292  const scalar time0,
293  const scalar time1,
294  const bool twoD,
295  const scalar angleOfWedge
296 ) const
297 {
298  scalar mInj = mass_*(fractionOfInjection(time1)-fractionOfInjection(time0));
299 
300  // correct mass if calculation is 2D
301  if (twoD)
302  {
303  mInj *= 0.5*angleOfWedge/mathematicalConstant::pi;
304  }
305 
306  return mInj;
307 }
308 
309 Foam::scalar Foam::multiHoleInjector::mass() const
310 {
311  return mass_;
312 }
313 
315 {
316  return X_;
317 }
318 
320 {
321  return TProfile_;
322 }
323 
324 Foam::scalar Foam::multiHoleInjector::T(const scalar time) const
325 {
326  return getTableValue(TProfile_, time);
327 }
328 
329 Foam::scalar Foam::multiHoleInjector::tsoi() const
330 {
331  return massFlowRateProfile_[0][0];
332 }
333 
334 Foam::scalar Foam::multiHoleInjector::teoi() const
335 {
336  return massFlowRateProfile_[massFlowRateProfile_.size()-1][0];
337 }
338 
340 (
341  const scalar time
342 ) const
343 {
344  return getTableValue(massFlowRateProfile_, time);
345 }
346 
348 (
349  const scalar time
350 ) const
351 {
352  return getTableValue(injectionPressureProfile_, time);
353 }
354 
356 (
357  const scalar time
358 ) const
359 {
360  return getTableValue(velocityProfile_, time);
361 }
362 
364 {
365  return CdProfile_;
366 }
367 
368 Foam::scalar Foam::multiHoleInjector::Cd
369 (
370  const scalar time
371 ) const
372 {
373  return Cd_;
374 }
375 
376 Foam::scalar Foam::multiHoleInjector::fractionOfInjection(const scalar time) const
377 {
378  return integrateTable(massFlowRateProfile_, time)/mass_;
379 }
380 
382 (
383  const scalar t
384 ) const
385 {
386  return mass_*fractionOfInjection(t);
387 }
388 
389 
391 (
392  const liquidMixture& fuel,
393  const scalar referencePressure
394 )
395 {
396 
397  scalar A = nHoles_*0.25*mathematicalConstant::pi*pow(d_, 2.0);
398 
399  forAll(velocityProfile_, i)
400  {
401  scalar time = velocityProfile_[i][0];
402  scalar rho = fuel.rho(referencePressure, T(time), X_);
403  scalar v = massFlowRateProfile_[i][1]/(Cd_*rho*A);
404  velocityProfile_[i][1] = v;
405  injectionPressureProfile_[i][1] = referencePressure + 0.5*rho*v*v;
406  }
407 }
408 
410 {
411  return tangentialInjectionVector1_[n];
412 }
413 
415 {
416  return tangentialInjectionVector2_[n];
417 }
418 
419 // ************************ vim: set sw=4 sts=4 et: ************************ //