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