src/lagrangian/dieselSpray/spray/sprayFunctions.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 1991-2009 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 "spray.H"
  28 #include "mathematicalConstants.H"
  29
  30 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
  31
  32 namespace Foam
  33 {
  34
  35 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
  36
  37 scalar spray::injectedMass(const scalar t) const
  38 {
  39     scalar sum = 0.0;
  40
  41     forAll (injectors_, i)
  42     {
  43         sum += injectors_[i].properties()->injectedMass(t);
  44     }
  45
  46     return sum;
  47 }
  48
  49
  50 scalar spray::totalMassToInject() const
  51 {
  52     scalar sum = 0.0;
  53
  54     forAll (injectors_, i)
  55     {
  56         sum += injectors_[i].properties()->mass();
  57     }
  58
  59     return sum;
  60 }
  61
  62
  63 scalar spray::injectedEnthalpy
  64 (
  65     const scalar time
  66 ) const
  67 {
  68     scalar sum = 0.0;
  69     label Nf = fuels_->components().size();
  70
  71     forAll (injectors_, i)
  72     {
  73         scalar T = injectors_[i].properties()->T(time);
  74         scalarField X(injectors_[i].properties()->X());
  75         scalar pi = 1.0e+5;
  76         scalar hl = fuels_->hl(pi, T, X);
  77         scalar Wl = fuels_->W(X);
  78         scalar hg = 0.0;
  79
  80         for(label j=0; j<Nf; j++)
  81         {
  82             label k = liquidToGasIndex_[j];
  83             hg += gasProperties()[k].H(T)*gasProperties()[k].W()*X[j]/Wl;
  84         }
  85
  86         sum += injectors_[i].properties()->injectedMass(time)*(hg-hl);
  87     }
  88
  89     return sum;
  90 }
  91
  92
  93 scalar spray::liquidMass() const
  94 {
  95     scalar sum = 0.0;
  96
  97     for
  98     (
  99         spray::const_iterator elmnt = begin();
 100         elmnt != end();
 101         ++elmnt
 102     )
 103     {
 104         sum += elmnt().m();
 105     }
 106
 107     if (twoD())
 108     {
 109         sum *= 2.0*mathematicalConstant::pi/angleOfWedge();
 110     }
 111
 112     reduce(sum, sumOp<scalar>());
 113
 114     return sum;
 115 }
 116
 117
 118 scalar spray::liquidEnthalpy() const
 119 {
 120     scalar sum = 0.0;
 121     label Nf = fuels().components().size();
 122
 123     for
 124     (
 125         spray::const_iterator elmnt = begin();
 126         elmnt != end();
 127         ++elmnt
 128     )
 129     {
 130         scalar T = elmnt().T();
 131         scalar pc = p()[elmnt().cell()];
 132         scalar hlat = fuels().hl(pc, T, elmnt().X());
 133         scalar hg = 0.0;
 134         scalar Wl = fuels().W(elmnt().X());
 135
 136         for(label j=0; j<Nf; j++)
 137         {
 138             label k = liquidToGasIndex_[j];
 139
 140             hg +=
 141                 gasProperties()[k].H(T)*gasProperties()[k].W()*elmnt().X()[j]
 142                /Wl;
 143         }
 144
 145         scalar h = hg - hlat;
 146         sum += elmnt().m()*h;
 147     }
 148
 149     if (twoD())
 150     {
 151         sum *= 2.0*mathematicalConstant::pi/angleOfWedge();
 152     }
 153
 154     reduce(sum, sumOp<scalar>());
 155
 156     return sum;
 157 }
 158
 159
 160 scalar spray::liquidTotalEnthalpy() const
 161 {
 162     scalar sum = 0.0;
 163     label Nf = fuels().components().size();
 164
 165     for
 166     (
 167         spray::const_iterator elmnt = begin();
 168         elmnt != end();
 169         ++elmnt
 170     )
 171     {
 172         label celli = elmnt().cell();
 173         scalar T = elmnt().T();
 174         scalar pc = p()[celli];
 175         scalar rho = fuels().rho(pc, T, elmnt().X());
 176         scalar hlat = fuels().hl(pc, T, elmnt().X());
 177         scalar hg = 0.0;
 178         scalar Wl = fuels().W(elmnt().X());
 179
 180         for(label j=0; j<Nf; j++)
 181         {
 182             label k = liquidToGasIndex_[j];
 183             hg +=
 184                 gasProperties()[k].H(T)*gasProperties()[k].W()*elmnt().X()[j]
 185                /Wl;
 186         }
 187
 188         scalar psat = fuels().pv(pc, T, elmnt().X());
 189
 190         scalar h = hg - hlat + (pc - psat)/rho;
 191         sum += elmnt().m()*h;
 192     }
 193
 194     if (twoD())
 195     {
 196         sum *= 2.0*mathematicalConstant::pi/angleOfWedge();
 197     }
 198
 199     reduce(sum, sumOp<scalar>());
 200
 201     return sum;
 202 }
 203
 204
 205 scalar spray::liquidKineticEnergy() const
 206 {
 207     scalar sum = 0.0;
 208     for
 209     (
 210         spray::const_iterator elmnt = begin();
 211         elmnt != end();
 212         ++elmnt
 213     )
 214     {
 215         scalar ke = pow(mag(elmnt().U()), 2.0);
 216         sum += elmnt().m()*ke;
 217     }
 218
 219     if (twoD())
 220     {
 221         sum *= 2.0*mathematicalConstant::pi/angleOfWedge();
 222     }
 223
 224     reduce(sum, sumOp<scalar>());
 225
 226     return 0.5*sum;
 227
 228 }
 229
 230
 231 scalar spray::injectedLiquidKineticEnergy() const
 232 {
 233     return injectedLiquidKE_;
 234 }
 235
 236
 237 scalar spray::liquidPenetration(const scalar prc) const
 238 {
 239     return liquidPenetration(0, prc);
 240 }
 241
 242
 243 scalar spray::liquidPenetration
 244 (
 245     const label nozzlei,
 246     const scalar prc
 247 ) const
 248 {
 249
 250     label nHoles = injectors_[nozzlei].properties()->nHoles();
 251     vector ip(vector::zero);
 252     if (nHoles > 1)
 253     {
 254         for(label i=0;i<nHoles;i++)
 255         {
 256             ip += injectors_[nozzlei].properties()->position(i);
 257         }
 258         ip /= nHoles;
 259     }
 260     else
 261     {
 262         ip = injectors_[nozzlei].properties()->position(0);
 263     }
 264
 265 //    vector ip = injectors_[nozzlei].properties()->position();
 266     scalar d = 0.0;
 267     scalar mTot = 0.0;
 268
 269     label Np = size();
 270
 271     // arrays containing the parcels mass and
 272     // distance from injector in ascending order
 273     scalarField m(Np);
 274     scalarField dist(Np);
 275     label n = 0;
 276
 277     if (Np > 1)
 278     {
 279         // NN.
 280         // first arrange the parcels in ascending order
 281         // the first parcel is closest to injector
 282         // and the last one is most far away.
 283         spray::const_iterator first = begin();
 284         m[n] = first().m();
 285         dist[n] = mag(first().position() - ip);
 286
 287         mTot += m[n];
 288
 289         for
 290         (
 291             spray::const_iterator elmnt = ++first;
 292             elmnt != end();
 293             ++elmnt
 294         )
 295         {
 296             scalar de = mag(elmnt().position() - ip);
 297             scalar me = elmnt().m();
 298             mTot += me;
 299
 300             n++;
 301
 302             label i = 0;
 303             bool found = false;
 304
 305             // insert the parcel in the correct place
 306             // and move the others
 307             while ( ( i < n-1 ) && ( !found ) )
 308             {
 309                 if (de < dist[i])
 310                 {
 311                     found = true;
 312                     for(label j=n; j>i; j--)
 313                     {
 314                         m[j]    = m[j-1];
 315                         dist[j] = dist[j-1];
 316                     }
 317                     m[i]    = me;
 318                     dist[i] = de;
 319                 }
 320                 i++;
 321             }
 322
 323             if (!found)
 324             {
 325                 m[n]    = me;
 326                 dist[n] = de;
 327             }
 328         }
 329     }
 330
 331     reduce(mTot, sumOp<scalar>());
 332
 333     if (Np > 1)
 334     {
 335         scalar mLimit = prc*mTot;
 336         scalar mOff = (1.0 - prc)*mTot;
 337
 338         // 'prc' is large enough that the parcel most far
 339         // away will be used, no need to loop...
 340         if (mLimit > mTot - m[Np-1])
 341         {
 342             d = dist[Np-1];
 343         }
 344         else
 345         {
 346             scalar mOffSum = 0.0;
 347             label i = Np;
 348
 349             while ((mOffSum < mOff) && (i>0))
 350             {
 351                 i--;
 352                 mOffSum += m[i];
 353             }
 354             d = dist[i];
 355         }
 356
 357     }
 358     else
 359     {
 360         if (Np > 0)
 361         {
 362             spray::const_iterator elmnt = begin();
 363             d = mag(elmnt().position() - ip);
 364         }
 365     }
 366
 367     reduce(d, maxOp<scalar>());
 368
 369     return d;
 370 }
 371
 372
 373 scalar spray::smd() const
 374 {
 375     scalar numerator = 0.0, denominator = VSMALL;
 376
 377     for
 378     (
 379         spray::const_iterator elmnt = begin();
 380         elmnt != end();
 381         ++elmnt
 382     )
 383     {
 384         label celli = elmnt().cell();
 385         scalar Pc = p()[celli];
 386         scalar T = elmnt().T();
 387         scalar rho = fuels_->rho(Pc, T, elmnt().X());
 388
 389         scalar tmp = elmnt().N(rho)*pow(elmnt().d(), 2.0);
 390         numerator += tmp*elmnt().d();
 391         denominator += tmp;
 392     }
 393
 394     reduce(numerator, sumOp<scalar>());
 395     reduce(denominator, sumOp<scalar>());
 396
 397     return numerator/denominator;
 398 }
 399
 400
 401 scalar spray::maxD() const
 402 {
 403     scalar maxD = 0.0;
 404
 405     for
 406     (
 407         spray::const_iterator elmnt = begin();
 408         elmnt != end();
 409         ++elmnt
 410     )
 411     {
 412         maxD = max(maxD, elmnt().d());
 413     }
 414
 415     reduce(maxD, maxOp<scalar>());
 416
 417     return maxD;
 418 }
 419
 420
 421 void spray::calculateAmbientPressure()
 422 {
 423     ambientPressure_ = p_.average().value();
 424 }
 425
 426
 427 void spray::calculateAmbientTemperature()
 428 {
 429     ambientTemperature_ = T_.average().value();
 430 }
 431
 432
 433 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 434
 435 } // End namespace Foam
 436
 437 // ************************************************************************* //