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