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