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) 1991-2008 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 the
  13     Free Software Foundation; either version 2 of the License, or (at your
  14     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, write to the Free Software Foundation,
  23     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
  24
  25 \*---------------------------------------------------------------------------*/
  26
  27 #include "unitInjector.H"
  28 #include "addToRunTimeSelectionTable.H"
  29 #include "Random.H"
  30 #include "mathematicalConstants.H"
  31
  32 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
  33 namespace Foam
  34 {
  35
  36 defineTypeNameAndDebug(unitInjector, 0);
  37
  38 addToRunTimeSelectionTable
  39 (
  40     injectorType,
  41     unitInjector,
  42     dictionary
  43 );
  44 }
  45 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  46
  47 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  48
  49 // Construct from components
  50 Foam::unitInjector::unitInjector
  51 (
  52     const Foam::Time& t,
  53     const Foam::dictionary& dict
  54 )
  55 :
  56     injectorType(t, dict),
  57     propsDict_(dict.subDict(typeName + "Props")),
  58     position_(propsDict_.lookup("position")),
  59     direction_(propsDict_.lookup("direction")),
  60     d_(readScalar(propsDict_.lookup("diameter"))),
  61     Cd_(readScalar(propsDict_.lookup("Cd"))),
  62     mass_(readScalar(propsDict_.lookup("mass"))),
  63     T_(readScalar(propsDict_.lookup("temperature"))),
  64     nParcels_(readLabel(propsDict_.lookup("nParcels"))),
  65     X_(propsDict_.lookup("X")),
  66     massFlowRateProfile_(propsDict_.lookup("massFlowRateProfile")),
  67     velocityProfile_(massFlowRateProfile_),
  68     injectionPressureProfile_(massFlowRateProfile_),
  69     CdProfile_(massFlowRateProfile_),
  70     TProfile_(massFlowRateProfile_),
  71     averageParcelMass_(mass_/nParcels_),
  72     pressureIndependentVelocity_(true)
  73 {
  74     // convert CA to real time
  75     forAll(massFlowRateProfile_, i)
  76     {
  77         massFlowRateProfile_[i][0] = t.userTimeToTime(massFlowRateProfile_[i][0]);
  78         velocityProfile_[i][0] = massFlowRateProfile_[i][0];
  79         injectionPressureProfile_[i][0] = massFlowRateProfile_[i][0];
  80     }
  81
  82     scalar integratedMFR = integrateTable(massFlowRateProfile_);
  83
  84     forAll(massFlowRateProfile_, i)
  85     {
  86         // correct the massFlowRateProfile to match the injected mass
  87         massFlowRateProfile_[i][1] *= mass_/integratedMFR;
  88
  89         TProfile_[i][0] = massFlowRateProfile_[i][0];
  90         TProfile_[i][1] = T_;
  91
  92         CdProfile_[i][0] = massFlowRateProfile_[i][0];
  93         CdProfile_[i][1] = Cd_;
  94     }
  95
  96     // Normalize the direction vector
  97     direction_ /= mag(direction_);
  98
  99     setTangentialVectors();
 100
 101     // check molar fractions
 102     scalar Xsum = 0.0;
 103     forAll(X_, i)
 104     {
 105         Xsum += X_[i];
 106     }
 107
 108     if (mag(Xsum - 1.0) > SMALL)
 109     {
 110         Info << "Warning!!!\n unitInjector::unitInjector(const time& t, Istream& is)"
 111             << "X does not add up to 1.0, correcting molar fractions."
 112             << endl;
 113         forAll(X_, i)
 114         {
 115             X_[i] /= Xsum;
 116         }
 117     }
 118
 119 }
 120
 121 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 122
 123 Foam::unitInjector::~unitInjector()
 124 {}
 125
 126
 127 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 128
 129 void Foam::unitInjector::setTangentialVectors()
 130 {
 131     Random rndGen(label(0));
 132     scalar magV = 0.0;
 133     vector tangent;
 134
 135     while (magV < SMALL)
 136     {
 137         vector testThis = rndGen.vector01();
 138
 139         tangent = testThis - (testThis & direction_)*direction_;
 140         magV = mag(tangent);
 141     }
 142
 143     tangentialInjectionVector1_ = tangent/magV;
 144     tangentialInjectionVector2_ = direction_ ^ tangentialInjectionVector1_;
 145
 146 }
 147
 148
 149 Foam::label Foam::unitInjector::nParcelsToInject
 150 (
 151     const scalar time0,
 152     const scalar time1
 153 ) const
 154 {
 155
 156     scalar mInj = mass_*(fractionOfInjection(time1)-fractionOfInjection(time0));
 157     label nParcels = label(mInj/averageParcelMass_ + 0.49);
 158
 159     return nParcels;
 160 }
 161
 162 const Foam::vector Foam::unitInjector::position() const
 163 {
 164     return position_;
 165 }
 166
 167 Foam::vector Foam::unitInjector::position
 168 (
 169     const scalar time,
 170     const bool twoD,
 171     const scalar angleOfWedge,
 172     const vector& axisOfSymmetry,
 173     const vector& axisOfWedge,
 174     const vector& axisOfWedgeNormal,
 175     Random& rndGen
 176 ) const
 177 {
 178     if (twoD)
 179     {
 180         scalar is = position_ & axisOfSymmetry;
 181         scalar magInj = mag(position_ - is*axisOfSymmetry);
 182
 183         vector halfWedge =
 184             axisOfWedge*cos(0.5*angleOfWedge)
 185           + axisOfWedgeNormal*sin(0.5*angleOfWedge);
 186         halfWedge /= mag(halfWedge);
 187
 188         return (is*axisOfSymmetry + magInj*halfWedge);
 189     }
 190     else
 191     {
 192         // otherwise, disc injection
 193         scalar iRadius = d_*rndGen.scalar01();
 194         scalar iAngle = 2.0*mathematicalConstant::pi*rndGen.scalar01();
 195
 196         return
 197         (
 198             position_
 199           + iRadius
 200           * (
 201               tangentialInjectionVector1_*cos(iAngle)
 202             + tangentialInjectionVector2_*sin(iAngle)
 203           )
 204         );
 205
 206     }
 207
 208     return position_;
 209 }
 210
 211 Foam::scalar Foam::unitInjector::d() const
 212 {
 213     return d_;
 214 }
 215
 216 const Foam::vector& Foam::unitInjector::direction() const
 217 {
 218     return direction_;
 219 }
 220
 221 Foam::scalar Foam::unitInjector::mass
 222 (
 223     const scalar time0,
 224     const scalar time1,
 225     const bool twoD,
 226     const scalar angleOfWedge
 227 ) const
 228 {
 229     scalar mInj = mass_*(fractionOfInjection(time1)-fractionOfInjection(time0));
 230
 231     // correct mass if calculation is 2D
 232     if (twoD)
 233     {
 234         mInj *= 0.5*angleOfWedge/mathematicalConstant::pi;
 235     }
 236
 237     return mInj;
 238 }
 239
 240 Foam::scalar Foam::unitInjector::mass() const
 241 {
 242     return mass_;
 243 }
 244
 245 const Foam::scalarField& Foam::unitInjector::X() const
 246 {
 247     return X_;
 248 }
 249
 250 Foam::List<Foam::unitInjector::pair> Foam::unitInjector::T() const
 251 {
 252     return TProfile_;
 253 }
 254
 255 Foam::scalar Foam::unitInjector::T(const scalar time) const
 256 {
 257     return T_;
 258 }
 259
 260 Foam::scalar Foam::unitInjector::tsoi() const
 261 {
 262     return massFlowRateProfile_[0][0];
 263 }
 264
 265 Foam::scalar Foam::unitInjector::teoi() const
 266 {
 267     return massFlowRateProfile_[massFlowRateProfile_.size()-1][0];
 268 }
 269
 270 Foam::scalar Foam::unitInjector::massFlowRate
 271 (
 272     const scalar time
 273 ) const
 274 {
 275     return getTableValue(massFlowRateProfile_, time);
 276 }
 277
 278 Foam::scalar Foam::unitInjector::injectionPressure
 279 (
 280     const scalar time
 281 ) const
 282 {
 283     return getTableValue(injectionPressureProfile_, time);
 284 }
 285
 286 Foam::scalar Foam::unitInjector::velocity
 287 (
 288     const scalar time
 289 ) const
 290 {
 291     return getTableValue(velocityProfile_, time);
 292 }
 293
 294 Foam::List<Foam::unitInjector::pair> Foam::unitInjector::CdProfile() const
 295 {
 296     return CdProfile_;
 297 }
 298
 299 Foam::scalar Foam::unitInjector::Cd
 300 (
 301     const scalar time
 302 ) const
 303 {
 304     return Cd_;
 305 }
 306
 307 Foam::scalar Foam::unitInjector::fractionOfInjection(const scalar time) const
 308 {
 309     return integrateTable(massFlowRateProfile_, time)/mass_;
 310 }
 311
 312 Foam::scalar Foam::unitInjector::injectedMass
 313 (
 314     const scalar t
 315 ) const
 316 {
 317     return mass_*fractionOfInjection(t);
 318 }
 319
 320
 321 void Foam::unitInjector::correctProfiles
 322 (
 323     const liquidMixture& fuel,
 324     const scalar referencePressure
 325 )
 326 {
 327
 328     scalar A = 0.25*mathematicalConstant::pi*pow(d_, 2.0);
 329     scalar pDummy = 1.0e+5;
 330
 331     scalar rho = fuel.rho(pDummy, T_, X_);
 332
 333     forAll(velocityProfile_, i)
 334     {
 335         scalar v = massFlowRateProfile_[i][1]/(Cd_*rho*A);
 336         velocityProfile_[i][1] = v;
 337         injectionPressureProfile_[i][1] = referencePressure + 0.5*rho*v*v;
 338     }
 339 }
 340
 341 // ************************************************************************* //