FreeFOAM The Cross-Platform CFD Toolkit
ThermoCloudI_.H
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) 1991-2010 OpenCFD Ltd.
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 
27 
28 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
29 
30 template<class ParcelType>
31 inline const typename ParcelType::constantProperties&
33 {
34  return constProps_;
35 }
36 
37 
38 template<class ParcelType>
39 inline const Foam::basicThermo&
41 {
42  return carrierThermo_;
43 }
44 
45 
46 template<class ParcelType>
47 inline Foam::basicThermo&
49 {
50  return carrierThermo_;
51 }
52 
53 
54 template<class ParcelType>
57 {
58  return heatTransferModel_;
59 }
60 
61 
62 template<class ParcelType>
65 {
66  return TIntegrator_;
67 }
68 
69 
70 template<class ParcelType>
72 {
73  return radiation_;
74 }
75 
76 
77 template<class ParcelType>
80 {
81  return hsTrans_;
82 }
83 
84 
85 template<class ParcelType>
88 {
90  (
92  (
93  IOobject
94  (
95  this->name() + "Sh",
96  this->db().time().timeName(),
97  this->mesh(),
98  IOobject::NO_READ,
99  IOobject::AUTO_WRITE,
100  false
101  ),
102  hsTrans_/(this->mesh().V()*this->db().time().deltaT())
103  )
104  );
105 
106  return tSh;
107 }
108 
109 
110 template<class ParcelType>
113 {
115  (
116  new volScalarField
117  (
118  IOobject
119  (
120  this->name() + "radiation::Ep",
121  this->db().time().timeName(),
122  this->db(),
123  IOobject::NO_READ,
124  IOobject::NO_WRITE,
125  false
126  ),
127  this->mesh(),
129  )
130  );
131 
132  // Need to check if coupled as field is created on-the-fly
133  if (radiation_ && this->coupled())
134  {
135  scalarField& Ep = tEp().internalField();
136  const scalarField& V = this->mesh().V();
137  const scalar epsilon = constProps_.epsilon0();
138 
139  forAllConstIter(typename ThermoCloud<ParcelType>, *this, iter)
140  {
141  const ParcelType& p = iter();
142  const label cellI = p.cell();
143  Ep[cellI] += p.nParticle()*p.areaP()*pow4(p.T());
144  }
145 
146  Ep *= epsilon*radiation::sigmaSB.value()/V;
147  }
148 
149  return tEp;
150 }
151 
152 
153 template<class ParcelType>
156 {
158  (
159  new volScalarField
160  (
161  IOobject
162  (
163  this->name() + "radiation::ap",
164  this->db().time().timeName(),
165  this->db(),
166  IOobject::NO_READ,
167  IOobject::NO_WRITE,
168  false
169  ),
170  this->mesh(),
171  dimensionedScalar("zero", dimless/dimLength, 0.0)
172  )
173  );
174 
175  // Need to check if coupled as field is created on-the-fly
176  if (radiation_ && this->coupled())
177  {
178  scalarField& ap = tap().internalField();
179  const scalarField& V = this->mesh().V();
180  const scalar epsilon = constProps_.epsilon0();
181 
182  forAllConstIter(typename ThermoCloud<ParcelType>, *this, iter)
183  {
184  const ParcelType& p = iter();
185  const label cellI = p.cell();
186  ap[cellI] += p.nParticle()*p.areaP();
187  }
188 
189  ap *= epsilon/V;
190  }
191 
192  return tap;
193 }
194 
195 
196 template<class ParcelType>
199 {
200  tmp<volScalarField> tsigmap
201  (
202  new volScalarField
203  (
204  IOobject
205  (
206  this->name() + "radiation::sigmap",
207  this->db().time().timeName(),
208  this->db(),
209  IOobject::NO_READ,
210  IOobject::NO_WRITE,
211  false
212  ),
213  this->mesh(),
214  dimensionedScalar("zero", dimless/dimLength, 0.0)
215  )
216  );
217 
218  // Need to check if coupled as field is created on-the-fly
219  if (radiation_ && this->coupled())
220  {
221  scalarField& sigmap = tsigmap().internalField();
222 
223  const scalarField& V = this->mesh().V();
224  const scalar epsilon = constProps_.epsilon0();
225  const scalar f = constProps_.f0();
226 
227  forAllConstIter(typename ThermoCloud<ParcelType>, *this, iter)
228  {
229  const ParcelType& p = iter();
230  const label cellI = p.cell();
231  sigmap[cellI] += p.nParticle()*p.areaP();
232  }
233 
234  sigmap *= (1.0 - f)*(1.0 - epsilon)/V;
235  }
236 
237  return tsigmap;
238 }
239 
240 
241 // ************************ vim: set sw=4 sts=4 et: ************************ //