src/lagrangian/dieselSpray/spray/sprayFunctions.C

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