src/lagrangian/dieselSpray/injector/unitInjector/unitInjector.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 2004-2011 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 "unitInjector.H"
  27 #include "addToRunTimeSelectionTable.H"
  28 #include "mathematicalConstants.H"
  29
  30 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  31
  32 namespace Foam
  33 {
  34     defineTypeNameAndDebug(unitInjector, 0);
  35
  36     addToRunTimeSelectionTable
  37     (
  38         injectorType,
  39         unitInjector,
  40         dictionary
  41     );
  42 }
  43
  44
  45 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  46
  47 Foam::unitInjector::unitInjector
  48 (
  49     const Foam::Time& t,
  50     const Foam::dictionary& dict
  51 )
  52 :
  53     injectorType(t, dict),
  54     propsDict_(dict.subDict(typeName + "Props")),
  55     position_(propsDict_.lookup("position")),
  56     direction_(propsDict_.lookup("direction")),
  57     d_(readScalar(propsDict_.lookup("diameter"))),
  58     Cd_(readScalar(propsDict_.lookup("Cd"))),
  59     mass_(readScalar(propsDict_.lookup("mass"))),
  60     nParcels_(readLabel(propsDict_.lookup("nParcels"))),
  61     X_(propsDict_.lookup("X")),
  62     massFlowRateProfile_(propsDict_.lookup("massFlowRateProfile")),
  63     velocityProfile_(massFlowRateProfile_),
  64     injectionPressureProfile_(massFlowRateProfile_),
  65     CdProfile_(massFlowRateProfile_),
  66     TProfile_(propsDict_.lookup("temperatureProfile")),
  67     averageParcelMass_(mass_/nParcels_),
  68     pressureIndependentVelocity_(true)
  69 {
  70
  71     // check if time entries for soi and eoi match
  72     if (mag(massFlowRateProfile_[0][0]-TProfile_[0][0]) > SMALL)
  73     {
  74         FatalErrorIn
  75         (
  76             "unitInjector::unitInjector(const time& t, const dictionary dict)"
  77         )   << "start-times do not match for TemperatureProfile and "
  78             << " massFlowRateProfile." << nl << exit (FatalError);
  79     }
  80
  81     if
  82     (
  83         mag(massFlowRateProfile_.last()[0]
  84       - TProfile_.last()[0])
  85       > SMALL
  86     )
  87     {
  88         FatalErrorIn
  89         (
  90             "unitInjector::unitInjector(const time& t, const dictionary dict)"
  91         )   << "end-times do not match for TemperatureProfile and "
  92             << "massFlowRateProfile." << nl << exit(FatalError);
  93     }
  94
  95     // convert CA to real time
  96     forAll(massFlowRateProfile_, i)
  97     {
  98         massFlowRateProfile_[i][0] =
  99             t.userTimeToTime(massFlowRateProfile_[i][0]);
 100         velocityProfile_[i][0] = massFlowRateProfile_[i][0];
 101         injectionPressureProfile_[i][0] = massFlowRateProfile_[i][0];
 102     }
 103
 104     forAll(TProfile_, i)
 105     {
 106         TProfile_[i][0] = t.userTimeToTime(TProfile_[i][0]);
 107     }
 108
 109     scalar integratedMFR = integrateTable(massFlowRateProfile_);
 110
 111     forAll(massFlowRateProfile_, i)
 112     {
 113         // correct the massFlowRateProfile to match the injected mass
 114         massFlowRateProfile_[i][1] *= mass_/integratedMFR;
 115
 116         CdProfile_[i][0] = massFlowRateProfile_[i][0];
 117         CdProfile_[i][1] = Cd_;
 118     }
 119
 120     // Normalize the direction vector
 121     direction_ /= mag(direction_);
 122
 123     setTangentialVectors();
 124
 125     // check molar fractions
 126     scalar Xsum = 0.0;
 127     forAll(X_, i)
 128     {
 129         Xsum += X_[i];
 130     }
 131
 132     if (mag(Xsum - 1.0) > SMALL)
 133     {
 134         WarningIn("unitInjector::unitInjector(const time& t, Istream& is)")
 135             << "X does not sum to 1.0, correcting molar fractions."
 136             << nl << endl;
 137         forAll(X_, i)
 138         {
 139             X_[i] /= Xsum;
 140         }
 141     }
 142
 143 }
 144
 145
 146 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 147
 148 Foam::unitInjector::~unitInjector()
 149 {}
 150
 151
 152 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 153
 154 void Foam::unitInjector::setTangentialVectors()
 155 {
 156     cachedRandom rndGen(label(0), -1);
 157     scalar magV = 0.0;
 158     vector tangent;
 159
 160     while (magV < SMALL)
 161     {
 162         vector testThis = rndGen.sample01<vector>();
 163
 164         tangent = testThis - (testThis & direction_)*direction_;
 165         magV = mag(tangent);
 166     }
 167
 168     tangentialInjectionVector1_ = tangent/magV;
 169     tangentialInjectionVector2_ = direction_ ^ tangentialInjectionVector1_;
 170 }
 171
 172
 173 Foam::label Foam::unitInjector::nParcelsToInject
 174 (
 175     const scalar time0,
 176     const scalar time1
 177 ) const
 178 {
 179     scalar mInj =
 180         mass_*(fractionOfInjection(time1) - fractionOfInjection(time0));
 181     label nParcels = label(mInj/averageParcelMass_ + 0.49);
 182     return nParcels;
 183 }
 184
 185
 186 const Foam::vector Foam::unitInjector::position(const label n) const
 187 {
 188     return position_;
 189 }
 190
 191
 192 Foam::vector Foam::unitInjector::position
 193 (
 194     const label n,
 195     const scalar time,
 196     const bool twoD,
 197     const scalar angleOfWedge,
 198     const vector& axisOfSymmetry,
 199     const vector& axisOfWedge,
 200     const vector& axisOfWedgeNormal,
 201     cachedRandom& rndGen
 202 ) const
 203 {
 204     if (twoD)
 205     {
 206         scalar is = position_ & axisOfSymmetry;
 207         scalar magInj = mag(position_ - is*axisOfSymmetry);
 208
 209         vector halfWedge =
 210             axisOfWedge*cos(0.5*angleOfWedge)
 211           + axisOfWedgeNormal*sin(0.5*angleOfWedge);
 212         halfWedge /= mag(halfWedge);
 213
 214         return (is*axisOfSymmetry + magInj*halfWedge);
 215     }
 216     else
 217     {
 218         // otherwise, disc injection
 219         scalar iRadius = d_*rndGen.sample01<scalar>();
 220         scalar iAngle = constant::mathematical::twoPi*rndGen.sample01<scalar>();
 221
 222         return
 223         (
 224             position_
 225           + iRadius
 226           * (
 227               tangentialInjectionVector1_*cos(iAngle)
 228             + tangentialInjectionVector2_*sin(iAngle)
 229           )
 230         );
 231     }
 232
 233     return position_;
 234 }
 235
 236
 237 Foam::label Foam::unitInjector::nHoles() const
 238 {
 239     return 1;
 240 }
 241
 242
 243 Foam::scalar Foam::unitInjector::d() const
 244 {
 245     return d_;
 246 }
 247
 248
 249 const Foam::vector& Foam::unitInjector::direction
 250 (
 251     const label i,
 252     const scalar time
 253 ) const
 254 {
 255     return direction_;
 256 }
 257
 258
 259 Foam::scalar Foam::unitInjector::mass
 260 (
 261     const scalar time0,
 262     const scalar time1,
 263     const bool twoD,
 264     const scalar angleOfWedge
 265 ) const
 266 {
 267     scalar mInj =
 268         mass_*(fractionOfInjection(time1) - fractionOfInjection(time0));
 269
 270     // correct mass if calculation is 2D
 271     if (twoD)
 272     {
 273         mInj *= 0.5*angleOfWedge/constant::mathematical::pi;
 274     }
 275
 276     return mInj;
 277 }
 278
 279
 280 Foam::scalar Foam::unitInjector::mass() const
 281 {
 282     return mass_;
 283 }
 284
 285
 286 const Foam::scalarField& Foam::unitInjector::X() const
 287 {
 288     return X_;
 289 }
 290
 291
 292 Foam::List<Foam::unitInjector::pair> Foam::unitInjector::T() const
 293 {
 294     return TProfile_;
 295 }
 296
 297
 298 Foam::scalar Foam::unitInjector::T(const scalar time) const
 299 {
 300     return getTableValue(TProfile_, time);
 301 }
 302
 303
 304 Foam::scalar Foam::unitInjector::tsoi() const
 305 {
 306     return massFlowRateProfile_.first()[0];
 307 }
 308
 309
 310 Foam::scalar Foam::unitInjector::teoi() const
 311 {
 312     return massFlowRateProfile_.last()[0];
 313 }
 314
 315
 316 Foam::scalar Foam::unitInjector::massFlowRate(const scalar time) const
 317 {
 318     return getTableValue(massFlowRateProfile_, time);
 319 }
 320
 321
 322 Foam::scalar Foam::unitInjector::injectionPressure(const scalar time) const
 323 {
 324     return getTableValue(injectionPressureProfile_, time);
 325 }
 326
 327
 328 Foam::scalar Foam::unitInjector::velocity(const scalar time) const
 329 {
 330     return getTableValue(velocityProfile_, time);
 331 }
 332
 333
 334 Foam::List<Foam::unitInjector::pair> Foam::unitInjector::CdProfile() const
 335 {
 336     return CdProfile_;
 337 }
 338
 339
 340 Foam::scalar Foam::unitInjector::Cd(const scalar time) const
 341 {
 342     return Cd_;
 343 }
 344
 345
 346 Foam::scalar Foam::unitInjector::fractionOfInjection(const scalar time) const
 347 {
 348     return integrateTable(massFlowRateProfile_, time)/mass_;
 349 }
 350
 351
 352 Foam::scalar Foam::unitInjector::injectedMass(const scalar t) const
 353 {
 354     return mass_*fractionOfInjection(t);
 355 }
 356
 357
 358 void Foam::unitInjector::correctProfiles
 359 (
 360     const liquidMixtureProperties& fuel,
 361     const scalar referencePressure
 362 )
 363 {
 364     scalar A = 0.25*constant::mathematical::pi*sqr(d_);
 365     scalar pDummy = 1.0e+5;
 366
 367     forAll(velocityProfile_, i)
 368     {
 369         scalar time = velocityProfile_[i][0];
 370         scalar rho = fuel.rho(pDummy, T(time), X_);
 371         scalar v = massFlowRateProfile_[i][1]/(Cd_*rho*A);
 372         velocityProfile_[i][1] = v;
 373         injectionPressureProfile_[i][1] = referencePressure + 0.5*rho*v*v;
 374     }
 375 }
 376
 377
 378 Foam::vector Foam::unitInjector::tan1(const label) const
 379 {
 380     return tangentialInjectionVector1_;
 381 }
 382
 383
 384 Foam::vector Foam::unitInjector::tan2(const label) const
 385 {
 386     return tangentialInjectionVector2_;
 387 }
 388
 389
 390 // ************************************************************************* //