FreeFOAM The Cross-Platform CFD Toolkit
blackBodyEmission.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) 2008-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 
26 #include "blackBodyEmission.H"
28 
29 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
30 
33 (
34  IStringStream
35  (
36  "("
37  "( 1000 0.00032)"
38  "( 1100 0.00091)"
39  "( 1200 0.00213)"
40  "( 1300 0.00432)"
41  "( 1400 0.00779)"
42  "( 1500 0.01280)"
43  "( 1600 0.01972)"
44  "( 1700 0.02853)"
45  "( 1800 0.03934)"
46  "( 1900 0.05210)"
47  "( 2000 0.06672)"
48  "( 2100 0.08305)"
49  "( 2200 0.10088)"
50  "( 2300 0.12002)"
51  "( 2400 0.14025)"
52  "( 2500 0.16135)"
53  "( 2600 0.18311)"
54  "( 2700 0.20535)"
55  "( 2800 0.22788)"
56  "( 2900 0.25055)"
57  "( 3000 0.27322)"
58  "( 3100 0.29576)"
59  "( 3200 0.31809)"
60  "( 3300 0.34009)"
61  "( 3400 0.36172)"
62  "( 3500 0.38290)"
63  "( 3600 0.40359)"
64  "( 3700 0.42375)"
65  "( 3800 0.44336)"
66  "( 3900 0.46240)"
67  "( 4000 0.48085)"
68  "( 4100 0.49872)"
69  "( 4200 0.51599)"
70  "( 4300 0.53267)"
71  "( 4400 0.54877)"
72  "( 4500 0.56429)"
73  "( 4600 0.57925)"
74  "( 4700 0.59366)"
75  "( 4800 0.60753)"
76  "( 4900 0.62088)"
77  "( 5000 0.63372)"
78  "( 5100 0.64606)"
79  "( 5200 0.65794)"
80  "( 5300 0.66935)"
81  "( 5400 0.68033)"
82  "( 5500 0.69087)"
83  "( 5600 0.70101)"
84  "( 5700 0.71076)"
85  "( 5800 0.72012)"
86  "( 5900 0.72913)"
87  "( 6000 0.73778)"
88  "( 6100 0.74610)"
89  "( 6200 0.75410)"
90  "( 6300 0.76180)"
91  "( 6400 0.76920)"
92  "( 6500 0.77631)"
93  "( 6600 0.78316)"
94  "( 6700 0.78975)"
95  "( 6800 0.79609)"
96  "( 6900 0.80219)"
97  "( 7000 0.80807)"
98  "( 7100 0.81373)"
99  "( 7200 0.81918)"
100  "( 7300 0.82443)"
101  "( 7400 0.82949)"
102  "( 7500 0.83436)"
103  "( 7600 0.83906)"
104  "( 7700 0.84359)"
105  "( 7800 0.84796)"
106  "( 7900 0.85218)"
107  "( 8000 0.85625)"
108  "( 8100 0.86017)"
109  "( 8200 0.86396)"
110  "( 8300 0.86762)"
111  "( 8400 0.87115)"
112  "( 8500 0.87456)"
113  "( 8600 0.87786)"
114  "( 8700 0.88105)"
115  "( 8800 0.88413)"
116  "( 8900 0.88711)"
117  "( 9000 0.88999)"
118  "( 9100 0.89277)"
119  "( 9200 0.89547)"
120  "( 9300 0.89807)"
121  "( 9400 0.90060)"
122  "( 9500 0.90304)"
123  "( 9600 0.90541)"
124  "( 9700 0.90770)"
125  "( 9800 0.90992)"
126  "( 9900 0.91207)"
127  "(10000 0.91415)"
128  "(12000 0.94505)"
129  "(15000 0.96893)"
130  "(20000 0.98555)"
131  "(30000 0.99529)"
132  "(40000 0.99792)"
133  "(50000 0.99890)"
134  ")"
135  )()
136 );
137 
138 
139 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
140 
142 (
143  const label nLambda,
144  const volScalarField& T
145 )
146 :
147  table_
148  (
149  emissivePowerTable,
151  "blackBodyEmissivePower"
152  ),
153  C1_("C1", dimensionSet(1, 4, 3, 0, 0, 0, 0), 3.7419e-16),
154  C2_("C2", dimensionSet(0, 1, 0, 1, 0, 0, 0), 14.388e-6),
155  bLambda_(nLambda),
156  T_(T)
157 {
158  forAll(bLambda_, lambdaI)
159  {
160  bLambda_.set
161  (
162  lambdaI,
163  new volScalarField
164  (
165  IOobject
166  (
167  "bLambda_" + Foam::name(lambdaI) ,
168  T.mesh().time().timeName(),
169  T.mesh(),
170  IOobject::NO_READ,
171  IOobject::NO_WRITE
172  ),
174  )
175  );
176 
177  }
178 }
179 
180 
181 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
182 
184 {}
185 
186 
187 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
188 
189 Foam::scalar Foam::radiation::blackBodyEmission::fLambdaT
190 (
191  const scalar lambdaT
192 ) const
193 {
194  return table_(lambdaT*1.0e6);
195 }
196 
197 
200 (
201  const volScalarField& T,
202  const Vector2D<scalar>& band
203 ) const
204 {
206  (
207  new volScalarField
208  (
209  IOobject
210  (
211  "Eb",
212  T.mesh().time().timeName(),
213  T.mesh(),
216  ),
218  )
219  );
220 
221 
222  if (band == Vector2D<scalar>::one)
223  {
224  return Eb;
225  }
226  else
227  {
228  forAll(T, i)
229  {
230  scalar T1 = fLambdaT(band[1]*T[i]);
231  scalar T2 = fLambdaT(band[0]*T[i]);
232  dimensionedScalar fLambdaDelta
233  (
234  "fLambdaDelta",
235  dimless,
236  T1 - T2
237  );
238  Eb()[i] = Eb()[i]*fLambdaDelta.value();
239  }
240  return Eb;
241  }
242 }
243 
244 
246 (
247  const label lambdaI,
248  const Vector2D<scalar>& band
249 )
250 {
251  bLambda_[lambdaI] = EbDeltaLambdaT(T_, band);
252 }
253 
254 
255 // ************************ vim: set sw=4 sts=4 et: ************************ //