src/thermophysicalModels/reactionThermo/combustionThermo/mixtureThermos/hhuMixtureThermo/hhuMixtureThermo.C

   1 /*---------------------------------------------------------------------------*\
   2   =========                 |
   3   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
   4    \\    /   O peration     |
   5     \\  /    A nd           | Copyright (C) 1991-2010 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
  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 "hhuMixtureThermo.H"
  27 #include "fvMesh.H"
  28 #include "fixedValueFvPatchFields.H"
  29
  30 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
  31
  32 template<class MixtureType>
  33 Foam::hhuMixtureThermo<MixtureType>::hhuMixtureThermo(const fvMesh& mesh)
  34 :
  35     hhuCombustionThermo(mesh),
  36     MixtureType(*this, mesh)
  37 {
  38     scalarField& hCells = h_.internalField();
  39     scalarField& huCells = hu_.internalField();
  40     const scalarField& TCells = T_.internalField();
  41     const scalarField& TuCells = Tu_.internalField();
  42
  43     forAll(hCells, celli)
  44     {
  45         hCells[celli] = this->cellMixture(celli).H(TCells[celli]);
  46         huCells[celli] = this->cellReactants(celli).H(TuCells[celli]);
  47     }
  48
  49     forAll(h_.boundaryField(), patchi)
  50     {
  51         h_.boundaryField()[patchi] == h(T_.boundaryField()[patchi], patchi);
  52
  53         fvPatchScalarField& phu = hu_.boundaryField()[patchi];
  54         const fvPatchScalarField& pTu = Tu_.boundaryField()[patchi];
  55
  56         forAll(phu, facei)
  57         {
  58             phu[facei] = this->patchFaceReactants(patchi, facei).H(pTu[facei]);
  59         }
  60     }
  61
  62     hBoundaryCorrection(h_);
  63     huBoundaryCorrection(hu_);
  64
  65     calculate();
  66     psi_.oldTime();   // Switch on saving old time
  67 }
  68
  69
  70 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
  71
  72 template<class MixtureType>
  73 Foam::hhuMixtureThermo<MixtureType>::~hhuMixtureThermo()
  74 {}
  75
  76
  77 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
  78
  79 template<class MixtureType>
  80 void Foam::hhuMixtureThermo<MixtureType>::calculate()
  81 {
  82     const scalarField& hCells = h_.internalField();
  83     const scalarField& huCells = hu_.internalField();
  84     const scalarField& pCells = p_.internalField();
  85
  86     scalarField& TCells = T_.internalField();
  87     scalarField& TuCells = Tu_.internalField();
  88     scalarField& psiCells = psi_.internalField();
  89     scalarField& muCells = mu_.internalField();
  90     scalarField& alphaCells = alpha_.internalField();
  91
  92     forAll(TCells, celli)
  93     {
  94         const typename MixtureType::thermoType& mixture_ =
  95             this->cellMixture(celli);
  96
  97         TCells[celli] = mixture_.TH(hCells[celli], TCells[celli]);
  98         psiCells[celli] = mixture_.psi(pCells[celli], TCells[celli]);
  99
 100         muCells[celli] = mixture_.mu(TCells[celli]);
 101         alphaCells[celli] = mixture_.alpha(TCells[celli]);
 102
 103         TuCells[celli] =
 104             this->cellReactants(celli).TH(huCells[celli], TuCells[celli]);
 105     }
 106
 107     forAll(T_.boundaryField(), patchi)
 108     {
 109         fvPatchScalarField& pp = p_.boundaryField()[patchi];
 110         fvPatchScalarField& pT = T_.boundaryField()[patchi];
 111         fvPatchScalarField& pTu = Tu_.boundaryField()[patchi];
 112         fvPatchScalarField& ppsi = psi_.boundaryField()[patchi];
 113
 114         fvPatchScalarField& ph = h_.boundaryField()[patchi];
 115         fvPatchScalarField& phu = hu_.boundaryField()[patchi];
 116
 117         fvPatchScalarField& pmu_ = mu_.boundaryField()[patchi];
 118         fvPatchScalarField& palpha_ = alpha_.boundaryField()[patchi];
 119
 120         if (pT.fixesValue())
 121         {
 122             forAll(pT, facei)
 123             {
 124                 const typename MixtureType::thermoType& mixture_ =
 125                     this->patchFaceMixture(patchi, facei);
 126
 127                 ph[facei] = mixture_.H(pT[facei]);
 128
 129                 ppsi[facei] = mixture_.psi(pp[facei], pT[facei]);
 130                 pmu_[facei] = mixture_.mu(pT[facei]);
 131                 palpha_[facei] = mixture_.alpha(pT[facei]);
 132             }
 133         }
 134         else
 135         {
 136             forAll(pT, facei)
 137             {
 138                 const typename MixtureType::thermoType& mixture_ =
 139                     this->patchFaceMixture(patchi, facei);
 140
 141                 pT[facei] = mixture_.TH(ph[facei], pT[facei]);
 142
 143                 ppsi[facei] = mixture_.psi(pp[facei], pT[facei]);
 144                 pmu_[facei] = mixture_.mu(pT[facei]);
 145                 palpha_[facei] = mixture_.alpha(pT[facei]);
 146
 147                 pTu[facei] =
 148                     this->patchFaceReactants(patchi, facei)
 149                    .TH(phu[facei], pTu[facei]);
 150             }
 151         }
 152     }
 153 }
 154
 155
 156 template<class MixtureType>
 157 void Foam::hhuMixtureThermo<MixtureType>::correct()
 158 {
 159     if (debug)
 160     {
 161         Info<< "entering hhuMixtureThermo<MixtureType>::correct()" << endl;
 162     }
 163
 164     // force the saving of the old-time values
 165     psi_.oldTime();
 166
 167     calculate();
 168
 169     if (debug)
 170     {
 171         Info<< "exiting hhuMixtureThermo<MixtureType>::correct()" << endl;
 172     }
 173 }
 174
 175
 176 template<class MixtureType>
 177 Foam::tmp<Foam::volScalarField>
 178 Foam::hhuMixtureThermo<MixtureType>::hc() const
 179 {
 180     const fvMesh& mesh = T_.mesh();
 181
 182     tmp<volScalarField> thc
 183     (
 184         new volScalarField
 185         (
 186             IOobject
 187             (
 188                 "hc",
 189                 mesh.time().timeName(),
 190                 mesh,
 191                 IOobject::NO_READ,
 192                 IOobject::NO_WRITE
 193             ),
 194             mesh,
 195             h_.dimensions()
 196         )
 197     );
 198
 199     volScalarField& hcf = thc();
 200     scalarField& hcCells = hcf.internalField();
 201
 202     forAll(hcCells, celli)
 203     {
 204         hcCells[celli] = this->cellMixture(celli).Hc();
 205     }
 206
 207     forAll(hcf.boundaryField(), patchi)
 208     {
 209         scalarField& hcp = hcf.boundaryField()[patchi];
 210
 211         forAll(hcp, facei)
 212         {
 213             hcp[facei] = this->patchFaceMixture(patchi, facei).Hc();
 214         }
 215     }
 216
 217     return thc;
 218 }
 219
 220
 221 template<class MixtureType>
 222 Foam::tmp<Foam::scalarField> Foam::hhuMixtureThermo<MixtureType>::h
 223 (
 224     const scalarField& T,
 225     const labelList& cells
 226 ) const
 227 {
 228     tmp<scalarField> th(new scalarField(T.size()));
 229     scalarField& h = th();
 230
 231     forAll(T, celli)
 232     {
 233         h[celli] = this->cellMixture(cells[celli]).H(T[celli]);
 234     }
 235
 236     return th;
 237 }
 238
 239
 240 template<class MixtureType>
 241 Foam::tmp<Foam::scalarField> Foam::hhuMixtureThermo<MixtureType>::h
 242 (
 243     const scalarField& T,
 244     const label patchi
 245 ) const
 246 {
 247     tmp<scalarField> th(new scalarField(T.size()));
 248     scalarField& h = th();
 249
 250     forAll(T, facei)
 251     {
 252         h[facei] = this->patchFaceMixture(patchi, facei).H(T[facei]);
 253     }
 254
 255     return th;
 256 }
 257
 258
 259 template<class MixtureType>
 260 Foam::tmp<Foam::scalarField> Foam::hhuMixtureThermo<MixtureType>::Cp
 261 (
 262     const scalarField& T,
 263     const label patchi
 264 ) const
 265 {
 266     tmp<scalarField> tCp(new scalarField(T.size()));
 267
 268     scalarField& cp = tCp();
 269
 270     forAll(T, facei)
 271     {
 272         cp[facei] = this->patchFaceMixture(patchi, facei).Cp(T[facei]);
 273     }
 274
 275     return tCp;
 276 }
 277
 278
 279 template<class MixtureType>
 280 Foam::tmp<Foam::volScalarField>
 281 Foam::hhuMixtureThermo<MixtureType>::Cp() const
 282 {
 283     const fvMesh& mesh = T_.mesh();
 284
 285     tmp<volScalarField> tCp
 286     (
 287         new volScalarField
 288         (
 289             IOobject
 290             (
 291                 "Cp",
 292                 mesh.time().timeName(),
 293                 mesh,
 294                 IOobject::NO_READ,
 295                 IOobject::NO_WRITE
 296             ),
 297             mesh,
 298             dimensionSet(0, 2, -2, -1, 0)
 299         )
 300     );
 301
 302     volScalarField& cp = tCp();
 303     scalarField& cpCells = cp.internalField();
 304     const scalarField& TCells = T_.internalField();
 305
 306     forAll(TCells, celli)
 307     {
 308         cpCells[celli] = this->cellMixture(celli).Cp(TCells[celli]);
 309     }
 310
 311     forAll(T_.boundaryField(), patchi)
 312     {
 313         cp.boundaryField()[patchi] = Cp(T_.boundaryField()[patchi], patchi);
 314     }
 315
 316     return tCp;
 317 }
 318
 319
 320 template<class MixtureType>
 321 Foam::tmp<Foam::scalarField>
 322 Foam::hhuMixtureThermo<MixtureType>::hu
 323 (
 324     const scalarField& Tu,
 325     const labelList& cells
 326 ) const
 327 {
 328     tmp<scalarField> thu(new scalarField(Tu.size()));
 329     scalarField& hu = thu();
 330
 331     forAll(Tu, celli)
 332     {
 333         hu[celli] = this->cellReactants(cells[celli]).H(Tu[celli]);
 334     }
 335
 336     return thu;
 337 }
 338
 339
 340 template<class MixtureType>
 341 Foam::tmp<Foam::scalarField>
 342 Foam::hhuMixtureThermo<MixtureType>::hu
 343 (
 344     const scalarField& Tu,
 345     const label patchi
 346 ) const
 347 {
 348     tmp<scalarField> thu(new scalarField(Tu.size()));
 349     scalarField& hu = thu();
 350
 351     forAll(Tu, facei)
 352     {
 353         hu[facei] = this->patchFaceReactants(patchi, facei).H(Tu[facei]);
 354     }
 355
 356     return thu;
 357 }
 358
 359
 360 template<class MixtureType>
 361 Foam::tmp<Foam::volScalarField>
 362 Foam::hhuMixtureThermo<MixtureType>::Tb() const
 363 {
 364     tmp<volScalarField> tTb
 365     (
 366         new volScalarField
 367         (
 368             IOobject
 369             (
 370                 "Tb",
 371                 T_.time().timeName(),
 372                 T_.db(),
 373                 IOobject::NO_READ,
 374                 IOobject::NO_WRITE
 375             ),
 376             T_
 377         )
 378     );
 379
 380     volScalarField& Tb_ = tTb();
 381     scalarField& TbCells = Tb_.internalField();
 382     const scalarField& TCells = T_.internalField();
 383     const scalarField& hCells = h_.internalField();
 384
 385     forAll(TbCells, celli)
 386     {
 387         TbCells[celli] =
 388             this->cellProducts(celli).TH(hCells[celli], TCells[celli]);
 389     }
 390
 391     forAll(Tb_.boundaryField(), patchi)
 392     {
 393         fvPatchScalarField& pTb = Tb_.boundaryField()[patchi];
 394
 395         const fvPatchScalarField& ph = h_.boundaryField()[patchi];
 396         const fvPatchScalarField& pT = T_.boundaryField()[patchi];
 397
 398         forAll(pTb, facei)
 399         {
 400             pTb[facei] =
 401                 this->patchFaceProducts(patchi, facei)
 402                .TH(ph[facei], pT[facei]);
 403         }
 404     }
 405
 406     return tTb;
 407 }
 408
 409
 410 template<class MixtureType>
 411 Foam::tmp<Foam::volScalarField>
 412 Foam::hhuMixtureThermo<MixtureType>::psiu() const
 413 {
 414     tmp<volScalarField> tpsiu
 415     (
 416         new volScalarField
 417         (
 418             IOobject
 419             (
 420                 "psiu",
 421                 psi_.time().timeName(),
 422                 psi_.db(),
 423                 IOobject::NO_READ,
 424                 IOobject::NO_WRITE
 425             ),
 426             psi_.mesh(),
 427             psi_.dimensions()
 428         )
 429     );
 430
 431     volScalarField& psiu = tpsiu();
 432     scalarField& psiuCells = psiu.internalField();
 433     const scalarField& TuCells = Tu_.internalField();
 434     const scalarField& pCells = p_.internalField();
 435
 436     forAll(psiuCells, celli)
 437     {
 438         psiuCells[celli] =
 439             this->cellReactants(celli).psi(pCells[celli], TuCells[celli]);
 440     }
 441
 442     forAll(psiu.boundaryField(), patchi)
 443     {
 444         fvPatchScalarField& ppsiu = psiu.boundaryField()[patchi];
 445
 446         const fvPatchScalarField& pp = p_.boundaryField()[patchi];
 447         const fvPatchScalarField& pTu = Tu_.boundaryField()[patchi];
 448
 449         forAll(ppsiu, facei)
 450         {
 451             ppsiu[facei] =
 452                 this->
 453                 patchFaceReactants(patchi, facei).psi(pp[facei], pTu[facei]);
 454         }
 455     }
 456
 457     return tpsiu;
 458 }
 459
 460
 461 template<class MixtureType>
 462 Foam::tmp<Foam::volScalarField>
 463 Foam::hhuMixtureThermo<MixtureType>::psib() const
 464 {
 465     tmp<volScalarField> tpsib
 466     (
 467         new volScalarField
 468         (
 469             IOobject
 470             (
 471                 "psib",
 472                 psi_.time().timeName(),
 473                 psi_.db(),
 474                 IOobject::NO_READ,
 475                 IOobject::NO_WRITE
 476             ),
 477             psi_.mesh(),
 478             psi_.dimensions()
 479         )
 480     );
 481
 482     volScalarField& psib = tpsib();
 483     scalarField& psibCells = psib.internalField();
 484     volScalarField Tb_ = Tb();
 485     const scalarField& TbCells = Tb_.internalField();
 486     const scalarField& pCells = p_.internalField();
 487
 488     forAll(psibCells, celli)
 489     {
 490         psibCells[celli] =
 491             this->cellReactants(celli).psi(pCells[celli], TbCells[celli]);
 492     }
 493
 494     forAll(psib.boundaryField(), patchi)
 495     {
 496         fvPatchScalarField& ppsib = psib.boundaryField()[patchi];
 497
 498         const fvPatchScalarField& pp = p_.boundaryField()[patchi];
 499         const fvPatchScalarField& pTb = Tb_.boundaryField()[patchi];
 500
 501         forAll(ppsib, facei)
 502         {
 503             ppsib[facei] =
 504                 this->patchFaceReactants
 505                 (patchi, facei).psi(pp[facei], pTb[facei]);
 506         }
 507     }
 508
 509     return tpsib;
 510 }
 511
 512
 513 template<class MixtureType>
 514 Foam::tmp<Foam::volScalarField>
 515 Foam::hhuMixtureThermo<MixtureType>::muu() const
 516 {
 517     tmp<volScalarField> tmuu
 518     (
 519         new volScalarField
 520         (
 521             IOobject
 522             (
 523                 "muu",
 524                 T_.time().timeName(),
 525                 T_.db(),
 526                 IOobject::NO_READ,
 527                 IOobject::NO_WRITE
 528             ),
 529             T_.mesh(),
 530             dimensionSet(1, -1, -1, 0, 0)
 531         )
 532     );
 533
 534     volScalarField& muu_ = tmuu();
 535     scalarField& muuCells = muu_.internalField();
 536     const scalarField& TuCells = Tu_.internalField();
 537
 538     forAll(muuCells, celli)
 539     {
 540         muuCells[celli] = this->cellReactants(celli).mu(TuCells[celli]);
 541     }
 542
 543     forAll(muu_.boundaryField(), patchi)
 544     {
 545         fvPatchScalarField& pMuu = muu_.boundaryField()[patchi];
 546         const fvPatchScalarField& pTu = Tu_.boundaryField()[patchi];
 547
 548         forAll(pMuu, facei)
 549         {
 550             pMuu[facei] =
 551                 this->patchFaceReactants(patchi, facei).mu(pTu[facei]);
 552         }
 553     }
 554
 555     return tmuu;
 556 }
 557
 558
 559 template<class MixtureType>
 560 Foam::tmp<Foam::volScalarField>
 561 Foam::hhuMixtureThermo<MixtureType>::mub() const
 562 {
 563     tmp<volScalarField> tmub
 564     (
 565         new volScalarField
 566         (
 567             IOobject
 568             (
 569                 "mub",
 570                 T_.time().timeName(),
 571                 T_.db(),
 572                 IOobject::NO_READ,
 573                 IOobject::NO_WRITE
 574             ),
 575             T_.mesh(),
 576             dimensionSet(1, -1, -1, 0, 0)
 577         )
 578     );
 579
 580     volScalarField& mub_ = tmub();
 581     scalarField& mubCells = mub_.internalField();
 582     volScalarField Tb_ = Tb();
 583     const scalarField& TbCells = Tb_.internalField();
 584
 585     forAll(mubCells, celli)
 586     {
 587         mubCells[celli] = this->cellProducts(celli).mu(TbCells[celli]);
 588     }
 589
 590     forAll(mub_.boundaryField(), patchi)
 591     {
 592         fvPatchScalarField& pMub = mub_.boundaryField()[patchi];
 593         const fvPatchScalarField& pTb = Tb_.boundaryField()[patchi];
 594
 595         forAll(pMub, facei)
 596         {
 597             pMub[facei] =
 598                 this->patchFaceProducts(patchi, facei).mu(pTb[facei]);
 599         }
 600     }
 601
 602     return tmub;
 603 }
 604
 605
 606 template<class MixtureType>
 607 bool Foam::hhuMixtureThermo<MixtureType>::read()
 608 {
 609     if (hhuCombustionThermo::read())
 610     {
 611         MixtureType::read(*this);
 612         return true;
 613     }
 614     else
 615     {
 616         return false;
 617     }
 618 }
 619
 620
 621 // ************************************************************************* //