Go to the documentation of this file.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 <radiation/radiationConstants.H>
00027
00028
00029
00030 template<class ParcelType>
00031 inline const typename ParcelType::constantProperties&
00032 Foam::ThermoCloud<ParcelType>::constProps() const
00033 {
00034 return constProps_;
00035 }
00036
00037
00038 template<class ParcelType>
00039 inline const Foam::basicThermo&
00040 Foam::ThermoCloud<ParcelType>::carrierThermo() const
00041 {
00042 return carrierThermo_;
00043 }
00044
00045
00046 template<class ParcelType>
00047 inline Foam::basicThermo&
00048 Foam::ThermoCloud<ParcelType>::carrierThermo()
00049 {
00050 return carrierThermo_;
00051 }
00052
00053
00054 template<class ParcelType>
00055 inline const Foam::HeatTransferModel<Foam::ThermoCloud<ParcelType> >&
00056 Foam::ThermoCloud<ParcelType>::heatTransfer() const
00057 {
00058 return heatTransferModel_;
00059 }
00060
00061
00062 template<class ParcelType>
00063 inline const Foam::scalarIntegrationScheme&
00064 Foam::ThermoCloud<ParcelType>::TIntegrator() const
00065 {
00066 return TIntegrator_;
00067 }
00068
00069
00070 template<class ParcelType>
00071 inline bool Foam::ThermoCloud<ParcelType>::radiation() const
00072 {
00073 return radiation_;
00074 }
00075
00076
00077 template<class ParcelType>
00078 inline Foam::DimensionedField<Foam::scalar, Foam::volMesh>&
00079 Foam::ThermoCloud<ParcelType>::hsTrans()
00080 {
00081 return hsTrans_;
00082 }
00083
00084
00085 template<class ParcelType>
00086 inline Foam::tmp<Foam::DimensionedField<Foam::scalar, Foam::volMesh> >
00087 Foam::ThermoCloud<ParcelType>::Sh() const
00088 {
00089 tmp<DimensionedField<scalar, volMesh> > tSh
00090 (
00091 new DimensionedField<scalar, volMesh>
00092 (
00093 IOobject
00094 (
00095 this->name() + "Sh",
00096 this->db().time().timeName(),
00097 this->mesh(),
00098 IOobject::NO_READ,
00099 IOobject::AUTO_WRITE,
00100 false
00101 ),
00102 hsTrans_/(this->mesh().V()*this->db().time().deltaT())
00103 )
00104 );
00105
00106 return tSh;
00107 }
00108
00109
00110 template<class ParcelType>
00111 inline Foam::tmp<Foam::volScalarField>
00112 Foam::ThermoCloud<ParcelType>::Ep() const
00113 {
00114 tmp<volScalarField> tEp
00115 (
00116 new volScalarField
00117 (
00118 IOobject
00119 (
00120 this->name() + "radiation::Ep",
00121 this->db().time().timeName(),
00122 this->db(),
00123 IOobject::NO_READ,
00124 IOobject::NO_WRITE,
00125 false
00126 ),
00127 this->mesh(),
00128 dimensionedScalar("zero", dimMass/dimLength/pow3(dimTime), 0.0)
00129 )
00130 );
00131
00132
00133 if (radiation_ && this->coupled())
00134 {
00135 scalarField& Ep = tEp().internalField();
00136 const scalarField& V = this->mesh().V();
00137 const scalar epsilon = constProps_.epsilon0();
00138
00139 forAllConstIter(typename ThermoCloud<ParcelType>, *this, iter)
00140 {
00141 const ParcelType& p = iter();
00142 const label cellI = p.cell();
00143 Ep[cellI] += p.nParticle()*p.areaP()*pow4(p.T());
00144 }
00145
00146 Ep *= epsilon*radiation::sigmaSB.value()/V;
00147 }
00148
00149 return tEp;
00150 }
00151
00152
00153 template<class ParcelType>
00154 inline Foam::tmp<Foam::volScalarField>
00155 Foam::ThermoCloud<ParcelType>::ap() const
00156 {
00157 tmp<volScalarField> tap
00158 (
00159 new volScalarField
00160 (
00161 IOobject
00162 (
00163 this->name() + "radiation::ap",
00164 this->db().time().timeName(),
00165 this->db(),
00166 IOobject::NO_READ,
00167 IOobject::NO_WRITE,
00168 false
00169 ),
00170 this->mesh(),
00171 dimensionedScalar("zero", dimless/dimLength, 0.0)
00172 )
00173 );
00174
00175
00176 if (radiation_ && this->coupled())
00177 {
00178 scalarField& ap = tap().internalField();
00179 const scalarField& V = this->mesh().V();
00180 const scalar epsilon = constProps_.epsilon0();
00181
00182 forAllConstIter(typename ThermoCloud<ParcelType>, *this, iter)
00183 {
00184 const ParcelType& p = iter();
00185 const label cellI = p.cell();
00186 ap[cellI] += p.nParticle()*p.areaP();
00187 }
00188
00189 ap *= epsilon/V;
00190 }
00191
00192 return tap;
00193 }
00194
00195
00196 template<class ParcelType>
00197 inline Foam::tmp<Foam::volScalarField>
00198 Foam::ThermoCloud<ParcelType>::sigmap() const
00199 {
00200 tmp<volScalarField> tsigmap
00201 (
00202 new volScalarField
00203 (
00204 IOobject
00205 (
00206 this->name() + "radiation::sigmap",
00207 this->db().time().timeName(),
00208 this->db(),
00209 IOobject::NO_READ,
00210 IOobject::NO_WRITE,
00211 false
00212 ),
00213 this->mesh(),
00214 dimensionedScalar("zero", dimless/dimLength, 0.0)
00215 )
00216 );
00217
00218
00219 if (radiation_ && this->coupled())
00220 {
00221 scalarField& sigmap = tsigmap().internalField();
00222
00223 const scalarField& V = this->mesh().V();
00224 const scalar epsilon = constProps_.epsilon0();
00225 const scalar f = constProps_.f0();
00226
00227 forAllConstIter(typename ThermoCloud<ParcelType>, *this, iter)
00228 {
00229 const ParcelType& p = iter();
00230 const label cellI = p.cell();
00231 sigmap[cellI] += p.nParticle()*p.areaP();
00232 }
00233
00234 sigmap *= (1.0 - f)*(1.0 - epsilon)/V;
00235 }
00236
00237 return tsigmap;
00238 }
00239
00240
00241