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