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-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 "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     vector ip = injectors_[nozzlei].properties()->position();
 250     scalar d = 0.0;
 251     scalar mTot = 0.0;
 252
 253     label Np = size();
 254
 255     // arrays containing the parcels mass and
 256     // distance from injector in ascending order
 257     scalarField m(Np);
 258     scalarField dist(Np);
 259     label n = 0;
 260
 261     if (Np > 1)
 262     {
 263         // NN.
 264         // first arrange the parcels in ascending order
 265         // the first parcel is closest to injector
 266         // and the last one is most far away.
 267         spray::const_iterator first = begin();
 268         m[n] = first().m();
 269         dist[n] = mag(first().position() - ip);
 270
 271         mTot += m[n];
 272
 273         for
 274         (
 275             spray::const_iterator elmnt = ++first;
 276             elmnt != end();
 277             ++elmnt
 278         )
 279         {
 280             scalar de = mag(elmnt().position() - ip);
 281             scalar me = elmnt().m();
 282             mTot += me;
 283
 284             n++;
 285
 286             label i = 0;
 287             bool found = false;
 288
 289             // insert the parcel in the correct place
 290             // and move the others
 291             while ( ( i < n-1 ) && ( !found ) )
 292             {
 293                 if (de < dist[i])
 294                 {
 295                     found = true;
 296                     for(label j=n; j>i; j--)
 297                     {
 298                         m[j]    = m[j-1];
 299                         dist[j] = dist[j-1];
 300                     }
 301                     m[i]    = me;
 302                     dist[i] = de;
 303                 }
 304                 i++;
 305             }
 306
 307             if (!found)
 308             {
 309                 m[n]    = me;
 310                 dist[n] = de;
 311             }
 312         }
 313     }
 314
 315     reduce(mTot, sumOp<scalar>());
 316
 317     if (Np > 1)
 318     {
 319         scalar mLimit = prc*mTot;
 320         scalar mOff = (1.0 - prc)*mTot;
 321
 322         // 'prc' is large enough that the parcel most far
 323         // away will be used, no need to loop...
 324         if (mLimit > mTot - m[Np-1])
 325         {
 326             d = dist[Np-1];
 327         }
 328         else
 329         {
 330             scalar mOffSum = 0.0;
 331             label i = Np;
 332
 333             while ((mOffSum < mOff) && (i>0))
 334             {
 335                 i--;
 336                 mOffSum += m[i];
 337             }
 338             d = dist[i];
 339         }
 340
 341     }
 342     else
 343     {
 344         if (Np > 0)
 345         {
 346             spray::const_iterator elmnt = begin();
 347             d = mag(elmnt().position() - ip);
 348         }
 349     }
 350
 351     reduce(d, maxOp<scalar>());
 352
 353     return d;
 354 }
 355
 356
 357 scalar spray::smd() const
 358 {
 359     scalar numerator = 0.0, denominator = VSMALL;
 360
 361     for
 362     (
 363         spray::const_iterator elmnt = begin();
 364         elmnt != end();
 365         ++elmnt
 366     )
 367     {
 368         label celli = elmnt().cell();
 369         scalar Pc = p()[celli];
 370         scalar T = elmnt().T();
 371         scalar rho = fuels_->rho(Pc, T, elmnt().X());
 372
 373         scalar tmp = elmnt().N(rho)*pow(elmnt().d(), 2.0);
 374         numerator += tmp*elmnt().d();
 375         denominator += tmp;
 376     }
 377
 378     reduce(numerator, sumOp<scalar>());
 379     reduce(denominator, sumOp<scalar>());
 380
 381     return numerator/denominator;
 382 }
 383
 384
 385 scalar spray::maxD() const
 386 {
 387     scalar maxD = 0.0;
 388
 389     for
 390     (
 391         spray::const_iterator elmnt = begin();
 392         elmnt != end();
 393         ++elmnt
 394     )
 395     {
 396         maxD = max(maxD, elmnt().d());
 397     }
 398
 399     reduce(maxD, maxOp<scalar>());
 400
 401     return maxD;
 402 }
 403
 404
 405 void spray::calculateAmbientPressure()
 406 {
 407     ambientPressure_ = p_.average().value();
 408 }
 409
 410
 411 void spray::calculateAmbientTemperature()
 412 {
 413     ambientTemperature_ = T_.average().value();
 414 }
 415
 416
 417 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 418
 419 } // End namespace Foam
 420
 421 // ************************************************************************* //