P1.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 2011-2016 OpenFOAM Foundation
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
13  the Free Software Foundation, either version 3 of the License, or
14  (at your 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, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "P1.H"
27 #include "fvmLaplacian.H"
28 #include "fvmSup.H"
30 #include "scatterModel.H"
31 #include "constants.H"
33 
34 using namespace Foam::constant;
35 
36 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
37 
38 namespace Foam
39 {
40  namespace radiation
41  {
44  }
45 }
46 
47 
48 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
49 
50 Foam::radiation::P1::P1(const volScalarField& T)
51 :
52  radiationModel(typeName, T),
53  G_
54  (
55  IOobject
56  (
57  "G",
58  mesh_.time().timeName(),
59  mesh_,
60  IOobject::MUST_READ,
61  IOobject::AUTO_WRITE
62  ),
63  mesh_
64  ),
65  Qr_
66  (
67  IOobject
68  (
69  "Qr",
70  mesh_.time().timeName(),
71  mesh_,
72  IOobject::NO_READ,
73  IOobject::AUTO_WRITE
74  ),
75  mesh_,
77  ),
78  a_
79  (
80  IOobject
81  (
82  "a",
83  mesh_.time().timeName(),
84  mesh_,
85  IOobject::NO_READ,
86  IOobject::AUTO_WRITE
87  ),
88  mesh_,
90  ),
91  e_
92  (
93  IOobject
94  (
95  "e",
96  mesh_.time().timeName(),
97  mesh_,
98  IOobject::NO_READ,
99  IOobject::NO_WRITE
100  ),
101  mesh_,
103  ),
104  E_
105  (
106  IOobject
107  (
108  "E",
109  mesh_.time().timeName(),
110  mesh_,
111  IOobject::NO_READ,
112  IOobject::NO_WRITE
113  ),
114  mesh_,
116  )
117 {}
118 
119 
120 Foam::radiation::P1::P1(const dictionary& dict, const volScalarField& T)
121 :
122  radiationModel(typeName, dict, T),
123  G_
124  (
125  IOobject
126  (
127  "G",
128  mesh_.time().timeName(),
129  mesh_,
132  ),
133  mesh_
134  ),
135  Qr_
136  (
137  IOobject
138  (
139  "Qr",
140  mesh_.time().timeName(),
141  mesh_,
142  IOobject::NO_READ,
144  ),
145  mesh_,
147  ),
148  a_
149  (
150  IOobject
151  (
152  "a",
153  mesh_.time().timeName(),
154  mesh_,
155  IOobject::NO_READ,
157  ),
158  mesh_,
160  ),
161  e_
162  (
163  IOobject
164  (
165  "e",
166  mesh_.time().timeName(),
167  mesh_,
168  IOobject::NO_READ,
170  ),
171  mesh_,
173  ),
174  E_
175  (
176  IOobject
177  (
178  "E",
179  mesh_.time().timeName(),
180  mesh_,
181  IOobject::NO_READ,
183  ),
184  mesh_,
186  )
187 {}
188 
189 
190 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
191 
193 {}
194 
195 
196 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
197 
199 {
200  if (radiationModel::read())
201  {
202  // nothing to read
203 
204  return true;
205  }
206  else
207  {
208  return false;
209  }
210 }
211 
212 
214 {
215  a_ = absorptionEmission_->a();
216  e_ = absorptionEmission_->e();
217  E_ = absorptionEmission_->E();
218  const volScalarField sigmaEff(scatter_->sigmaEff());
219 
220  const dimensionedScalar a0 ("a0", a_.dimensions(), ROOTVSMALL);
221 
222  // Construct diffusion
223  const volScalarField gamma
224  (
225  IOobject
226  (
227  "gammaRad",
228  G_.mesh().time().timeName(),
229  G_.mesh(),
232  ),
233  1.0/(3.0*a_ + sigmaEff + a0)
234  );
235 
236  // Solve G transport equation
237  solve
238  (
239  fvm::laplacian(gamma, G_)
240  - fvm::Sp(a_, G_)
241  ==
242  - 4.0*(e_*physicoChemical::sigma*pow4(T_) ) - E_
243  );
244 
245  volScalarField::Boundary& QrBf = Qr_.boundaryFieldRef();
246 
247  // Calculate radiative heat flux on boundaries.
249  {
250  if (!G_.boundaryField()[patchi].coupled())
251  {
252  QrBf[patchi] =
253  -gamma.boundaryField()[patchi]
254  *G_.boundaryField()[patchi].snGrad();
255  }
256  }
257 }
258 
259 
261 {
262  return tmp<volScalarField>
263  (
264  new volScalarField
265  (
266  IOobject
267  (
268  "Rp",
269  mesh_.time().timeName(),
270  mesh_,
273  false
274  ),
276  )
277  );
278 }
279 
280 
283 {
285  G_();
287  absorptionEmission_->ECont()()();
289  absorptionEmission_->aCont()()();
290 
291  return a*G - E;
292 }
293 
294 
295 // ************************************************************************* //
Collection of constants.
dictionary dict
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:428
autoPtr< scatterModel > scatter_
Scatter model.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:137
tmp< fvMatrix< Type > > Sp(const DimensionedField< scalar, volMesh > &, const GeometricField< Type, fvPatchField, volMesh > &)
const dimensionedScalar G
Newtonian constant of gravitation.
const volScalarField & T_
Reference to the temperature field.
const dimensionedScalar sigma
Stefan-Boltzmann constant: default SI units: [W/m2/K4].
Calculate the matrix for the laplacian of the field.
bool read()
Read radiation properties dictionary.
Definition: P1.C:198
static word timeName(const scalar, const int precision=precision_)
Return time name of given scalar time.
Definition: Time.C:715
const polyBoundaryMesh & boundaryMesh() const
Return boundary mesh.
Definition: polyMesh.H:421
Macros for easy insertion into run-time selection tables.
GeometricField< scalar, fvPatchField, volMesh > volScalarField
Definition: volFieldsFwd.H:52
virtual tmp< volScalarField > Rp() const
Source term component (for power of T^4)
Definition: P1.C:260
virtual bool read()=0
Read radiationProperties dictionary.
const Boundary & boundaryField() const
Return const-reference to the boundary field.
SolverPerformance< Type > solve(fvMatrix< Type > &, const dictionary &)
Solve returning the solution statistics given convergence tolerance.
const Time & time() const
Return time.
Definition: IOobject.C:227
void calculate()
Solve radiation equation(s)
Definition: P1.C:213
const fvMesh & mesh_
Reference to the mesh database.
Top level model for radiation modelling.
word timeName
Definition: getTimeIndex.H:3
virtual ~P1()
Destructor.
Definition: P1.C:192
const dimensionSet & dimensions() const
Return dimensions.
defineTypeNameAndDebug(combustionModel, 0)
void T(FieldField< Field, Type > &f1, const FieldField< Field, Type > &f2)
dimensionedScalar pow3(const dimensionedScalar &ds)
label patchi
const dimensionSet dimless(0, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:47
Boundary & boundaryFieldRef()
Return a reference to the boundary field.
const dimensionSet dimLength(0, 1, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:50
tmp< fvMatrix< Type > > laplacian(const GeometricField< Type, fvPatchField, volMesh > &vf, const word &name)
Definition: fvmLaplacian.C:46
const Mesh & mesh() const
Return mesh.
virtual tmp< DimensionedField< scalar, volMesh > > Ru() const
Source term component (constant)
Definition: P1.C:282
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensionedScalar pow4(const dimensionedScalar &ds)
const dimensionedScalar a0
Bohr radius: default SI units: [m].
const dimensionSet dimTime(0, 0, 1, 0, 0, 0, 0)
Definition: dimensionSets.H:51
const dimensionSet dimMass(1, 0, 0, 0, 0, 0, 0)
Definition: dimensionSets.H:49
A class for managing temporary objects.
Definition: PtrList.H:54
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:91
#define addToRadiationRunTimeSelectionTables(model)
autoPtr< absorptionEmissionModel > absorptionEmission_
Absorption/emission model.
Calculate the matrix for implicit and explicit sources.
autoPtr< radiation::radiationModel > radiation(radiation::radiationModel::New(T))
Namespace for OpenFOAM.
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:243