FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

P1.C

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software: you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
00023 
00024 \*---------------------------------------------------------------------------*/
00025 
00026 #include "P1.H"
00027 #include <OpenFOAM/addToRunTimeSelectionTable.H>
00028 #include <finiteVolume/fvm.H>
00029 #include <finiteVolume/fvc.H>
00030 
00031 #include <radiation/absorptionEmissionModel.H>
00032 #include <radiation/scatterModel.H>
00033 #include <OpenFOAM/mathematicalConstants.H>
00034 #include <radiation/radiationConstants.H>
00035 
00036 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
00037 
00038 namespace Foam
00039 {
00040     namespace radiation
00041     {
00042         defineTypeNameAndDebug(P1, 0);
00043 
00044         addToRunTimeSelectionTable
00045         (
00046             radiationModel,
00047             P1,
00048             dictionary
00049         );
00050     }
00051 }
00052 
00053 
00054 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00055 
00056 Foam::radiation::P1::P1(const volScalarField& T)
00057 :
00058     radiationModel(typeName, T),
00059     G_
00060     (
00061         IOobject
00062         (
00063             "G",
00064             mesh_.time().timeName(),
00065             mesh_,
00066             IOobject::MUST_READ,
00067             IOobject::AUTO_WRITE
00068         ),
00069         mesh_
00070     ),
00071     Qr_
00072     (
00073         IOobject
00074         (
00075             "Qr",
00076             mesh_.time().timeName(),
00077             mesh_,
00078             IOobject::NO_READ,
00079             IOobject::AUTO_WRITE
00080         ),
00081         mesh_,
00082         dimensionedScalar("Qr", dimMass/pow3(dimTime), 0.0)
00083     ),
00084     a_
00085     (
00086         IOobject
00087         (
00088             "a",
00089             mesh_.time().timeName(),
00090             mesh_,
00091             IOobject::NO_READ,
00092             IOobject::NO_WRITE
00093         ),
00094         mesh_,
00095         dimensionedScalar("a", dimless/dimLength, 0.0)
00096     ),
00097     e_
00098     (
00099         IOobject
00100         (
00101             "e",
00102             mesh_.time().timeName(),
00103             mesh_,
00104             IOobject::NO_READ,
00105             IOobject::NO_WRITE
00106         ),
00107         mesh_,
00108         dimensionedScalar("a", dimless/dimLength, 0.0)
00109     ),
00110     E_
00111     (
00112         IOobject
00113         (
00114             "E",
00115             mesh_.time().timeName(),
00116             mesh_,
00117             IOobject::NO_READ,
00118             IOobject::NO_WRITE
00119         ),
00120         mesh_,
00121         dimensionedScalar("E", dimMass/dimLength/pow3(dimTime), 0.0)
00122     )
00123 {}
00124 
00125 
00126 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
00127 
00128 Foam::radiation::P1::~P1()
00129 {}
00130 
00131 
00132 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00133 
00134 bool Foam::radiation::P1::read()
00135 {
00136     if (radiationModel::read())
00137     {
00138         // nothing to read
00139 
00140         return true;
00141     }
00142     else
00143     {
00144         return false;
00145     }
00146 }
00147 
00148 
00149 void Foam::radiation::P1::calculate()
00150 {
00151     a_ = absorptionEmission_->a();
00152     e_ = absorptionEmission_->e();
00153     E_ = absorptionEmission_->E();
00154     const volScalarField sigmaEff = scatter_->sigmaEff();
00155 
00156     // Construct diffusion
00157     const volScalarField gamma
00158     (
00159         IOobject
00160         (
00161             "gammaRad",
00162             G_.mesh().time().timeName(),
00163             G_.mesh(),
00164             IOobject::NO_READ,
00165             IOobject::NO_WRITE
00166         ),
00167         1.0/(3.0*a_ + sigmaEff)
00168     );
00169 
00170     // Solve G transport equation
00171     solve
00172     (
00173         fvm::laplacian(gamma, G_)
00174       - fvm::Sp(a_, G_)
00175      ==
00176       - 4.0*(e_*radiation::sigmaSB*pow4(T_) + E_)
00177     );
00178 
00179     // Calculate radiative heat flux on boundaries.
00180     forAll(mesh_.boundaryMesh(), patchI)
00181     {
00182         Qr_.boundaryField()[patchI] =
00183             -gamma.boundaryField()[patchI]*G_.boundaryField()[patchI].snGrad();
00184     }
00185 }
00186 
00187 
00188 Foam::tmp<Foam::volScalarField> Foam::radiation::P1::Rp() const
00189 {
00190     return tmp<volScalarField>
00191     (
00192         new volScalarField
00193         (
00194             IOobject
00195             (
00196                 "Rp",
00197                 mesh_.time().timeName(),
00198                 mesh_,
00199                 IOobject::NO_READ,
00200                 IOobject::NO_WRITE,
00201                 false
00202             ),
00203             4.0*absorptionEmission_->eCont()*radiation::sigmaSB
00204         )
00205     );
00206 }
00207 
00208 
00209 Foam::tmp<Foam::DimensionedField<Foam::scalar, Foam::volMesh> >
00210 Foam::radiation::P1::Ru() const
00211 {
00212     const DimensionedField<scalar, volMesh>& G =
00213         G_.dimensionedInternalField();
00214     const DimensionedField<scalar, volMesh> E =
00215         absorptionEmission_->ECont()().dimensionedInternalField();
00216     const DimensionedField<scalar, volMesh> a =
00217         absorptionEmission_->aCont()().dimensionedInternalField();
00218 
00219     return  a*G - 4.0*E;
00220 }
00221 
00222 
00223 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines