00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
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
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
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
00127
00128 Foam::radiation::P1::~P1()
00129 {}
00130
00131
00132
00133
00134 bool Foam::radiation::P1::read()
00135 {
00136 if (radiationModel::read())
00137 {
00138
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
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
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
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