src/finiteVolume/cfdTools/general/porousMedia/porousZoneTemplates.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 "porousZone.H"
  28 #include "fvMesh.H"
  29
  30 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
  31
  32 template<class Type>
  33 void Foam::porousZone::modifyDdt(fvMatrix<Type>& m) const
  34 {
  35     if (porosity_ < 1)
  36     {
  37         const labelList& cells = mesh_.cellZones()[cellZoneID_];
  38
  39         forAll(cells, i)
  40         {
  41             m.diag()[cells[i]]   *= porosity_;
  42             m.source()[cells[i]] *= porosity_;
  43         }
  44     }
  45 }
  46
  47
  48 template<class RhoFieldType>
  49 void Foam::porousZone::addPowerLawResistance
  50 (
  51     scalarField& Udiag,
  52     const labelList& cells,
  53     const scalarField& V,
  54     const RhoFieldType& rho,
  55     const vectorField& U
  56 ) const
  57 {
  58     const scalar C0 = C0_;
  59     const scalar C1m1b2 = (C1_ - 1.0)/2.0;
  60
  61     forAll (cells, i)
  62     {
  63         Udiag[cells[i]] +=
  64             V[cells[i]]*rho[cells[i]]*C0*pow(magSqr(U[cells[i]]), C1m1b2);
  65     }
  66 }
  67
  68
  69 template<class RhoFieldType>
  70 void Foam::porousZone::addViscousInertialResistance
  71 (
  72     scalarField& Udiag,
  73     vectorField& Usource,
  74     const labelList& cells,
  75     const scalarField& V,
  76     const RhoFieldType& rho,
  77     const scalarField& mu,
  78     const vectorField& U
  79 ) const
  80 {
  81     const tensor& D = D_.value();
  82     const tensor& F = F_.value();
  83
  84     forAll (cells, i)
  85     {
  86         tensor dragCoeff = mu[cells[i]]*D + (rho[cells[i]]*mag(U[cells[i]]))*F;
  87         scalar isoDragCoeff = tr(dragCoeff);
  88
  89         Udiag[cells[i]] += V[cells[i]]*isoDragCoeff;
  90         Usource[cells[i]] -=
  91             V[cells[i]]*((dragCoeff - I*isoDragCoeff) & U[cells[i]]);
  92     }
  93 }
  94
  95
  96 template<class RhoFieldType>
  97 void Foam::porousZone::addPowerLawResistance
  98 (
  99     tensorField& AU,
 100     const labelList& cells,
 101     const RhoFieldType& rho,
 102     const vectorField& U
 103 ) const
 104 {
 105     const scalar C0 = C0_;
 106     const scalar C1m1b2 = (C1_ - 1.0)/2.0;
 107
 108     forAll (cells, i)
 109     {
 110         AU[cells[i]] = AU[cells[i]]
 111           + I*(rho[cells[i]]*C0*pow(magSqr(U[cells[i]]), C1m1b2));
 112     }
 113 }
 114
 115
 116 template<class RhoFieldType>
 117 void Foam::porousZone::addViscousInertialResistance
 118 (
 119     tensorField& AU,
 120     const labelList& cells,
 121     const RhoFieldType& rho,
 122     const scalarField& mu,
 123     const vectorField& U
 124 ) const
 125 {
 126     const tensor& D = D_.value();
 127     const tensor& F = F_.value();
 128
 129     forAll (cells, i)
 130     {
 131         AU[cells[i]] += mu[cells[i]]*D + (rho[cells[i]]*mag(U[cells[i]]))*F;
 132     }
 133 }
 134
 135
 136 // ************************************************************************* //