00001 /*---------------------------------------------------------------------------*\ 00002 ========= | 00003 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox 00004 \\ / O peration | 00005 \\ / A nd | Copyright (C) 2008-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 Class 00025 Foam::radiation::fvDOM 00026 00027 Description 00028 00029 Finite Volume Discrete Ordinates Method. Solves the RTE equation for n 00030 directions in a participating media, not including scatter. 00031 00032 Available absorption models: 00033 greyMeanAbsoprtionEmission 00034 wideBandAbsorptionEmission 00035 00036 i.e. dictionary 00037 fvDOMCoeffs 00038 { 00039 nPhi 1; // azimuthal angles in PI/2 on X-Y.(from Y to X) 00040 nTheta 2; // polar angles in PI (from Z to X-Y plane) 00041 convergence 1e-4; // convergence criteria for radiation iteration 00042 } 00043 00044 solverFreq 1; // Number of flow iterations per radiation iteration 00045 00046 The total number of solid angles is 4*nPhi*nTheta. 00047 00048 In 1D the direction of the rays is X (nPhi and nTheta are ignored) 00049 In 2D the direction of the rays is on X-Y plane (only nPhi is considered) 00050 In 3D (nPhi and nTheta are considered) 00051 00052 SourceFiles 00053 fvDOM.C 00054 00055 \*---------------------------------------------------------------------------*/ 00056 00057 #ifndef radiationModelfvDOM_H 00058 #define radiationModelfvDOM_H 00059 00060 #include <radiation/radiativeIntensityRay.H> 00061 #include <radiation/radiationModel.H> 00062 00063 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 00064 00065 namespace Foam 00066 { 00067 namespace radiation 00068 { 00069 00070 /*---------------------------------------------------------------------------*\ 00071 Class fvDOM Declaration 00072 \*---------------------------------------------------------------------------*/ 00073 00074 class fvDOM 00075 : 00076 public radiationModel 00077 { 00078 // Private data 00079 00080 //- Incident radiation [W/m2] 00081 volScalarField G_; 00082 00083 //- Total radiative heat flux [W/m2] 00084 volScalarField Qr_; 00085 00086 //- Emmited radiative heat flux [W/m2] 00087 volScalarField Qem_; 00088 00089 //- Incidet radiative heat flux [W/m2] 00090 volScalarField Qin_; 00091 00092 //- Total absorption coefficient [1/m] 00093 volScalarField a_; 00094 00095 //- Total emission coefficient [1/m] 00096 volScalarField e_; 00097 00098 //- Emission contribution [Kg/m/s^3] 00099 volScalarField E_; 00100 00101 //- Number of solid angles in theta 00102 label nTheta_; 00103 00104 //- Number of solid angles in phi 00105 label nPhi_ ; 00106 00107 //- Total number of rays (1 per direction) 00108 label nRay_; 00109 00110 //- Number of wavelength bands 00111 label nLambda_; 00112 00113 //- Wavelength total absorption coefficient [1/m] 00114 PtrList<volScalarField> aLambda_; 00115 00116 //- Black body 00117 blackBodyEmission blackBody_; 00118 00119 //- List of pointers to radiative intensity rays 00120 PtrList<radiativeIntensityRay> IRay_; 00121 00122 //- Convergence criterion 00123 scalar convergence_; 00124 00125 //- Maximum number of iterations 00126 scalar maxIter_; 00127 00128 00129 // Private member functions 00130 00131 //- Disallow default bitwise copy construct 00132 fvDOM(const fvDOM&); 00133 00134 //- Disallow default bitwise assignment 00135 void operator=(const fvDOM&); 00136 00137 //- Update nlack body emission 00138 void updateBlackBodyEmission(); 00139 00140 00141 public: 00142 00143 //- Runtime type information 00144 TypeName("fvDOM"); 00145 00146 00147 // Constructors 00148 00149 //- Construct from components 00150 fvDOM(const volScalarField& T); 00151 00152 00153 //- Destructor 00154 virtual ~fvDOM(); 00155 00156 00157 // Member functions 00158 00159 // Edit 00160 00161 //- Solve radiation equation(s) 00162 void calculate(); 00163 00164 //- Read radiation properties dictionary 00165 bool read(); 00166 00167 //- Update G and calculate total heat flux on boundary 00168 void updateG(); 00169 00170 //- Set the rayId and lambdaId from by decomposing an intensity 00171 // field name 00172 void setRayIdLambdaId 00173 ( 00174 const word& name, 00175 label& rayId, 00176 label& lambdaId 00177 ) const; 00178 00179 //- Source term component (for power of T^4) 00180 virtual tmp<volScalarField> Rp() const; 00181 00182 //- Source term component (constant) 00183 virtual tmp<DimensionedField<scalar, volMesh> > Ru() const; 00184 00185 00186 // Access 00187 00188 //- Ray intensity for rayI 00189 inline const radiativeIntensityRay& IRay(const label rayI) const; 00190 00191 //- Ray intensity for rayI and lambda bandwidth 00192 inline const volScalarField& IRayLambda 00193 ( 00194 const label rayI, 00195 const label lambdaI 00196 ) const; 00197 00198 //- Number of angles in theta 00199 inline label nTheta() const; 00200 00201 //- Number of angles in phi 00202 inline label nPhi() const; 00203 00204 //- Number of rays 00205 inline label nRay() const; 00206 00207 //- Number of wavelengths 00208 inline label nLambda() const; 00209 00210 //- Const access to total absorption coefficient 00211 inline const volScalarField& a() const; 00212 00213 //- Const access to wavelength total absorption coefficient 00214 inline const volScalarField& aLambda(const label lambdaI) const; 00215 00216 //- Const access to incident radiation field 00217 inline const volScalarField& G() const; 00218 00219 //- Const access to total radiative heat flux field 00220 inline const volScalarField& Qr() const; 00221 00222 //- Const access to incident radiative heat flux field 00223 inline const volScalarField& Qin() const; 00224 00225 //- Const access to emitted radiative heat flux field 00226 inline const volScalarField& Qem() const; 00227 00228 //- Const access to black body 00229 inline const blackBodyEmission& blackBody() const; 00230 }; 00231 00232 00233 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 00234 00235 #include "fvDOMI.H" 00236 00237 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 00238 00239 } // End namespace radiation 00240 } // End namespace Foam 00241 00242 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 00243 00244 #endif 00245 00246 // ************************ vim: set sw=4 sts=4 et: ************************ //