0
Research Papers: Radiative Heat Transfer

Linearization of Nongray Radiation Exchange: The Internal Fractional Function Reconsidered OPEN ACCESS

[+] Author and Article Information
John H. Lienhard V

Rohsenow Kendall Heat Transfer Laboratory,
Department of Mechanical Engineering,
Massachusetts Institute of Technology,
Cambridge, MA 02139-4307
e-mail: lienhard@mit.edu

1The internal radiation fractional function also arises in models for radiation extinction in stellar interiors [, Sec. 37] and other optically thick, nongray gases [], but those applications are unrelated to the present discussion.

Contributed by the Heat Transfer Division of ASME for publication in the JOURNAL OF HEAT TRANSFER. Manuscript received July 4, 2018; final manuscript received November 20, 2018; published online March 27, 2019. Assoc. Editor: Milind A. Jog.

J. Heat Transfer 141(5), 052701 (Mar 27, 2019) (9 pages) Paper No: HT-18-1432; doi: 10.1115/1.4042158 History: Received July 04, 2018; Revised November 20, 2018

The radiation fractional function is the fraction of black body radiation below a given value of λT. Edwards and others have distinguished between the traditional, or “external,” radiation fractional function and an “internal” radiation fractional function. The latter is used for linearization of net radiation from a nongray surface when the temperature of an effectively black environment is not far from the surface's temperature, without calculating a separate total absorptivity. This paper examines the analytical approximation involved in the internal fractional function, with results given in terms of the incomplete zeta function. A rigorous upper bound on the difference between the external and internal emissivity is obtained. Calculations using the internal emissivity are compared to exact calculations for several models and materials. A new approach to calculating the internal emissivity is developed, yielding vastly improved accuracy over a wide range of temperature differences. The internal fractional function should be used for evaluating radiation thermal resistances, in particular.

FIGURES IN THIS ARTICLE
<>

The use of the radiation fractional function is well known in computing total emissivity from spectral emissivity. Less well known is the use of the internal radiation fractional function for approximating net heat flux from one surface to another at a similar temperature. Edwards was a strong proponent of the latter approach, beginning in the 1960s [1], following work published a few years earlier by Czerny and Walther [2]. Edwards' work on radiative property measurements had well acquainted him with the failure of gray-body exchange approximations at even small temperature differences [3], leading him to seek a better approximation that retained the computational convenience of using only a single radiative property.

Edwards distinguished the traditional fractional function by calling it the external fractional function. The terms external and internal derive from Edwards' work in the aerospace industry, differentiating between radiation exchange on the exterior of a spacecraft, where temperature differences are large, and on the interior of a spacecraft, where surface temperatures are not far apart [4]. The internal radiation fractional function appears in several textbooks by Edwards and coworkers [57] spanning a period of decades. Most other textbooks have not taken it up. In part, the lack of adoption may reflect the limited range of temperature difference over which the internal emissivity provides an accurate estimate of the net heat flux.

The present paper reconsiders the approximation leading to the internal radiation fractional function, with the objective of understanding its accuracy and how that may be improved. The results are given in terms of the incomplete zeta function, and a rigorous upper bound on the difference between the external and internal emissivity is obtained. Calculations based upon the internal radiation fractional function are compared to exact calculations, and a procedure is established that sharply increases domain of accuracy of the internal emissivity. The internal emissivity should be used in computing radiation thermal resistance in multimode heat transfer and for other linearizations of nongray radiation exchange.

The spectral hemispherical absorptivity, α(λ, T1), of a surface at T1, is related to the total hemispherical2 emissivity, ε(T1), by integration over the Planck distribution of blackbody monochromatic emissive power, eλ,b(T1) Display Formula

(1)ε(T1)=1σT140α(λ,T1)eλ,b(T1)dλ
Display Formula
(2)eλ,b(T1)=πc1λ5[exp(c2/λT1)1]

where c1 and c2 are the first and second radiation constants. If this surface is irradiated by an effectively black environment at T2 (for example, if the surroundings are much larger than the surface), the total hemispherical absorptivity is Display Formula

(3)α(T1,T2)=1σT240α(λ,T1)eλ,b(T2)dλ

The net radiation leaving this surface is Display Formula

(4)qnet=σε(T1)T14σα(T1,T2)T24

For a gray surface, α(λ, T) is independent of wavelength and so ε(T1) = α(T1, T2). However, some older literature went farther, applying this idea to nongray surfaces for which T1 is near T2. Edwards noted that the fact that limT2T1α(T1,T2)=ε(T1) is meaningless in calculating heat exchange; and further, working with Nelson, he documented associated “nongray errors” in radiative property measurements that exceeded 60% [3].

The reason for such errors is not that the limit is incorrect, but that the heat flux is determined by the slope of qnet as a function of temperature difference, which changes if the surface is not truly gray. To see why in a simple way, the last term in Eq. (4) may be linearized Display Formula

(5)α(T1,T2)T24α(T1,T1)T14+ddT2(α(T1,T2)T24)|T1(T2T1)
Display Formula
(6)=ε(T1)T14+4T13[ε(T1)+T14dαdT2|T1](T2T1)

Thus Display Formula

(7)qnet4σT13[ε(T1)+T14dαdT2|T1](T1T2)

For a gray (or black) surface, /dT2 = 0, and we recover a common linearization: qnet4σεT13ΔT [11]. Thus, the purpose of introducing the internal fractional function is to account for the difference in slope when /dT2 ≠ 0.

The fraction of blackbody radiation between wavelengths of 0 and λ is Display Formula

(8)f(λT)=1σT40λeλ,bdλ
Display Formula
(9)=115π40c2/λTt3et1dt
Display Formula
(10)=190π4ζ(c2/λT,4)

ζ(X,s) is the incomplete zeta function (see derivation in Appendix A). Edwards [1,4,5] refers to f as the external radiation fraction. The distinct internal radiation fractional function, fi(λT), arises from a linearization that allows ε(Ts) and α(Ts,Te) to be combined, as described next.

Edwards [1] defines the internal hemispherical emissivity such that Display Formula

(11)εi(T1)limT2T1σε(T1)T14σα(T1,T2)T24σT14σT24

Because the numerator and denominator both go to zero in the limit, L'Hôpital's rule (or Taylor expansion) is required Display Formula

(12)εi(T1)=limT2T10α(λ,T1)T2eλ,b(T2)dλ4σT23
Display Formula
(13)=0α(λ,T1)[14σT13T2eλ,b(T2)|T1]dλ

Thus, when T2 is not too much different from T1Display Formula

(14)qnetεi(T1)4σT13(T1T2)

with Display Formula

(15)εi(T)=14σT30α(λ,T)eλ,bTdλ=01α(λ,T)dfi(λT)

where the internal radiation fractional function is Display Formula

(16)fi(λT)14σT30λeλ,bTdλ

L'Hôpital's rule simply compares the slopes of the numerator and denominator in the limit, so Eqs. (7) and (14) are identical statements, as a change of variables and integration by parts of Eq. (13) shows.

Equation (15) stands in contrast to the external total hemispherical emissivity, which from Eqs. (1) and (8) maybe written Display Formula

(17)ε(T)=1σT40α(λ,T)eλ,b(T)dλ=01α(λ,T)df(λT)

From these relationships, one can show that Display Formula

(18)fi(λT)f(λT)=F(X)=154π4X4eX1

where Xc2/λT. The functions f, fi, and F are plotted as a function of λT in Fig. 1.

The numerical implementation is described in Appendix B.

We now derive an upper bound on the difference of these two emissivities Display Formula

(19)εεi=01α(λ,T)df(λT)01α(λ,T)dfi(λT)
Display Formula
(20)=0α(λ,T)dFdXdX
Display Formula
(21)=0Xzα(λ,T)dFdXdX+Xzα(λ,T)dFdXdX

where X = Xz is the finite zero point of dF/dX, Xz = 3.92069 (see Appendix C). The function dF/dX is positive for X < Xz and negative for X > Xz, and 0<α(λ,T)1. Thus, we can form two bounds on the difference Display Formula

(22)εεi0XzdFdXdX=F(Xz)ifεεi>0,and
Display Formula
(23)εiεXzdFdXdX=F(Xz)ifεiε>0

A simple calculation shows that Display Formula

(24)F(Xz)=15π4Xz3eXz=0.18400

In other words Display Formula

(25)|εεi|0.18400

The value of Xz corresponds to Display Formula

(26)(λT)z=14387.773.92069=3699.70μmK

At 300 K, the wavelength λz = 12.2323 μm.

The following two extreme examples are designed to illustrate these bounds. First, suppose that we have a 300 K surface that is black for λ12.2323μm, but perfectly reflective on other wavelengths. Numerical integration of Eqs. (15) and (17) yields Display Formula

(27)ε=0.4177,εi=0.6017,andεiε=0.1840

If instead we have a 300 K surface that is black for 12.2323μmλ, but perfectly reflective on other wavelengths Display Formula

(28)ε=0.5823,εi=0.3983,andεεi=0.1840

In the first case, εi is 44% greater than ε; in the second case, εi is 32% less than ε.

For a gray surface, Eq. (20) gives ε − εi = 0, consistent with Eqs. (4) and (7). The variation of the total absorptivity for a 300 K surface irradiated by an effectively black source at T2 is shown for both example surfaces in Fig. 2. The total absorptivity varies significantly. The net heat exchange is discussed in detail below.

The emissivities obtained from these examples are more general than they may appear: the same values will be obtained for any surface for which α(λ) switches between zero and one at Xz. For example, from the preceding development, for the surface black on short wavelengths Display Formula

(29)ε=0λzα(λ)df(λTz)=XzdfdXdX=f(Xz)0
Display Formula
(30)=190π4ζ(Xz,4)=190π4ζ(3.92069,4)
Display Formula
(31)=0.41771

Equation (18) gives εi, and the results for the long wavelength example may be obtained by subtracting these values from unity.

The linearization itself produces error in calculation of the heat flux, apart from nongray error. To quantify this error, we may compare the linearized and exact values of qnet for a black surface. Their ratio is Display Formula

(32)4T13ΔTT14T24=4T13(T12+T22)(T1+T2)=4[1+(T2/T1)2](1+T2/T1)

Equation (32) is plotted in Fig. 3. For T2 = 320 K and T1 = 300 K, the linearized flux is 90.5% of the exact flux. The ratio is asymptotic to 13ΔT/2T1+5(ΔT)2/4T12 as ΔT → 0.

The linearization around T1 is considerably less accurate than the traditional approach (i.e., for radiation heat transfer coefficients) of linearizing around Tm = (T1 + T2)/2. For that case, the ratio is Display Formula

(33)4Tm3ΔTT14T24=4Tm3(2Tm2+ΔT2/2)(2Tm)=11+ΔT2/4Tm2

and has errors on the order of (ΔT/Tm)2/4 [11]. As seen in Fig. 3, this linearization is far more accurate. The reason for the greater accuracy is that linearization about Tm is equivalent to a second-order Runge–Kutta method, with third-order truncation error in heat flux, while linearization about T1 is simply a forward Euler method (see Appendix D).

From these considerations, we see that εi will also have a limited range of accuracy if evaluated at T1. Instead, εi should be computed at Tm, following the analysis in Appendix D: Display Formula

(34)εi(Tm)=14σTm30α(λ,T1)eλ,bT|Tmdλ=01α(λ,T1)dfi(λTm)
with the heat flux evaluated, to an accuracy of O(ΔT3), as Display Formula
(35)qnet4εi(Tm)·σTm3ΔT

We compare the heat flux predicted by the exact result (Eq. (4)), the gray approximation, εi(T1) (Eq. (14)), and εi(Tm) (Eq. (35)) in Figs. 4 and 5 for the same two example surfaces at T1 = 300 K. The charts show

  1. (1)The gray surface approximation, ε = α, has the wrong slope as T2T1, consistent with Eq. (7). This approximation has poor accuracy for even the smallest temperature differences. The ratio of gray slope to exact slope is ε/εi. Using the values from Eq. (27) or Eq. (28), that ratio is equal to 0.69 and 1.46 for Figs. 4 and 5, respectively. As expected, the gray approximation shows the greatest error for T2 > T1, where irradiation dominates net exchange.
  2. (2)The approximation using εi(T1) has the correct slope as T2T1, but is accurate for only a limited range of temperatures. For the surface black on short wavelengths, the magnitude of the linearized flux from Eq. (14) at T2 = 320 K is 13% less than the exact value, Eq. (4). This case corresponds to the seemly small temperature difference ΔT/T1 = −0.067.
  3. (3)The approximation using εi(Tm) remains accurate over a very broad range of temperature variation, to the degree that the curves for this equation and the exact equation are difficult to distinguish. The charts in Figs. 4 and 5 include values as high as ΔT/T1 = ±0.33.

Although Figs. 4 and 5 use 300 K as a reference point, Eqs. (29)(31) show that, for a surface at T1 = Tz that switches from black to reflective behavior at λz, both ε(Tz) and εi(Tz) are fixed. The values of εi(Tm)=εi(TzΔT/2) and α(Tz, Tz − ΔT) depend additionally only upon ΔT/Tz. Hence, the trends seen will be the same at any other temperature level. These comparisons are put into nondimensional form in Appendix F.

Edwards [1,4,5] and, later, Mills [7] were both very specific in recommending to evaluate εi at T1. Edwards et al. [6] did not specify the temperature, but without comment included an example worked using Tm. Likewise, Czerny and Walther [2] provided an example worked using Tm, again without comment. The present results show clearly that only Tm should be used when evaluating εi.

The previous example surfaces were contrived to maximize the differences between ε and εi. In this section, we examine the differences for representative nongray materials.

Alumina.

Experimental data for polycrystalline alumina (99.5% Al2O3, 6 mm thick, 1 μm roughness) have been reported by Teodorescu and Jones [12] for T =823 K and 2λ25μm. The experimental uncertainty is 3.5%. These data are available in tabular form on 1 μm intervals and as a chart (Fig. 6). These results are generally consistent with other reports on polycrystalline alumina.

In this case, we first integrate the data to obtain the spectral normal emissivity, εn. One additional data point is extrapolated, ε(1 μm, 823 K) ≈ 0.128, because fi at 2 μm still has the finite value 0.059; the extrapolation extends the range to fi = 0.00012. At the other end of this spectrum, f =0.987 at 25 μm, but no extrapolation was done. Thus, the value for total normal emissivity εn is thus likely to be low by ∼1%. Integration of the discrete data by a simple trapezoidal rule yields the following values, both at 823 K: Display Formula

(36)εn=0.506,εni=0.404,andεniεn=0.102

Our interest, ultimately, lies with estimating heat exchange, so we require the spectral hemispherical emissivity. Teodorescu and Jones [12] also measured the spectral directional emissivity. We have integrated those data, as described in Appendix E, to obtain the spectral hemispherical values from 1 to 25 μm.

We consider exchange with an effectively black environment at temperatures down to 625 K. At 625 K, f =0.972 at 25 μm so that the absorptivity calculation may be low by 2% or so. The results are shown in Fig. 7. Once again, Eq. (35) using εi(Tm) is in excellent agreement with the exact result, Eq. (4).

Drude/Hagens-Ruben Metal.

The Drude/Hagen–Rubens model for emissivity applies approximately to some metals in the infrared regime [5,10,13]. Baehr and Stephan [13] provide the following equation for the spectral, hemispherical emissivity under this model: Display Formula

(37)ε(λ,T)=48.70reλ{1+[31.62+6.849ln(reλ)]reλ166.78reλ+}
for re/λ<5×104Ω·cm/μm, where re(T) is the electrical resistivity in Ω cm.

We consider a surface at T1 = 373 K with re = 13.1 × 10−6 Ω ⋅ cm, a value representative of platinum. At this temperature, 0.01f0.99 for 3.88λ61.3μm and fi = 0.01 at λ = 3.36 μm. Over this range of wavelengths, ε(λ) varies from 0.086 down to 0.022 according to Eq. (37). We may consider heat exchange with a surface at T2 down to 275 K; at 275 K, f >0.99 for λ > 83.1 μm. The integrations are, therefore, done from 3 to 100 μm. The results are shown in Fig. 8.

Once again, Eq. (35) using εi(Tm) closely tracks the exact solution over a wide range (down to at least 300 K or ΔT/T1 = 20%). In this instance, nongray effects are present but less pronounced; the gray approximation, while not the best alternative for even small ΔT, shows less divergence than in previous examples. In addition, at 373 K Display Formula

(38)ε=0.0497,εi=0.0553,andεiε=0.0056

The internal and external emissivities are closer, as would be expected when nongray effects are smaller.

Spectrally Selective Materials.

We may also consider an idealized spectrally selective surface that makes a step transition at λc from a short-wavelength absorptivity, αsw, to a long-wavelength absorptivity, αlwDisplay Formula

(39)α(λ)={αswforλλcαlwforλ>λc

From Eqs. (10) and (17)Display Formula

(40)ε(T)=αswf(λcT)+αlw[1f(λcT)]
Display Formula
(41)=αsw+90π4Δαζ(Xc,4)

where Xc = c2/λcT and Δα = αlw − αsw. Similarly, from Eqs. (15) and (18)Display Formula

(42)εi(T)=αswfi(λcT)+αlw[1fi(λcT)]
Display Formula
(43)=αsw+Δα[90π4ζ(Xc,4)F(Xc)]

Further Display Formula

(44)ε(T)εi(T)=ΔαF(Xc)

and also Display Formula

(45)εi(Tm)=αsw+Δα[90π4ζ(Xc,m,4)F(Xc,m)]

where Xc,m = c2/λcTm. Finally, from Eq. (3)Display Formula

(46)α(T1,T2)=αsw+90π4Δαζ(Xc,2,4)
with Xc,2 = c2/λcT2. The previous results show that the impact of selectivity will be greatest when Xc and Xz are close.

As an example, we may consider the following crude approximation to the soft-anodized aluminum described in Ref. [5]: αsw = 0.10, αlw = 0.85, and λc = 7 μm (Fig. 9). If this selective solar reflector is at T1 = 360 K and exchanges radiation with sky at T2 = 290 K, we obtain the results in Table 1. At these temperatures, most of the radiant energy is on wavelengths above 7 μm, and nongrayness is not pronounced. Nonetheless, as in previous examples, the agreement between Eq. (4), 371.8 W/m2, and Eq. (35), 371.0 W/m2, is excellent, whereas Eq. (14), 462.1 W/m2, performs quite poorly.

Equations (39)(46) are used further in Appendix F.

One rationale for introducing εi is simply to avoid computing the total hemispherical emissivity and absorptivity separately when the spectral absorptivity has a significant wavelength dependence. The same labor-saving rationale was used in some older literature to justify a gray body approximation, ε(T1) ≈ α(T1, T2), when the temperatures are close. When the internal fractional function was first developed, personal computers (and pocket calculators, spreadsheets, symbolic integration, etc.) did not exist. Routine analyses required either tedious hand calculations or cumbersome trips to the mainframe computer center (to type up Hollerith cards and then wait for the SysOp to run the job!). So, one attraction of the method was clearly to reduce labor.

While numerical integration is very easy today, εi(Tm) can still save some time in doing basic calculations. The greater challenge may be to obtain complete spectral, hemispherical data spanning the necessary range of wavelengths. If, instead, electromagnetic theory is used, working from the complex refractive index and calculating spectral directional emissivity, much coding has to be undertaken, and the single integration saved by εi will not be at all significant.

However, a second very important consideration is that a linearized expression for radiation heat transfer allows the use of a radiation thermal resistance, perhaps in parallel to a natural convection resistance, in network analyses of heat transfer. For an area A, such a resistance has the form Display Formula

(47)Rrad=14εσT3A

Rrad is routinely evaluated using a gray body approximation. For small temperature differences, using the internal emissivity, εi(T1), in Eq. (47) will provide a more accurate value than a gray body approximation. The present work has shown εi(Tm) to be far more robust, with excellent accuracy for |ΔT/T1| of 20–30%. Hence, the radiation thermal resistance is best calculated as Display Formula

(48)Rrad=14εi(Tm)σTm3A

Despite modern computing power, the need for an accurate radiation resistance remains as important today as it was 60 years ago. Thus, we conclude that the internal fractional function has an ongoing value when surfaces do not politely behave as “gray bodies.”

The internal radiation fractional function has been in the literature since the late 1950s, although not widely adopted in textbooks covering thermal radiation. The most visible recommendations have been from Edwards and coworkers, generally advising the use of εi(T1) to linearize nongray exchange for enclosures with modest temperature differences. In this study, the following new findings are reported:

  1. (1)The maximum absolute difference between ε(T1) and εi(T1) is 0.18400, Eq. (25).
  2. (2)The internal emissivity should be evaluated at the mean temperature, Tm, not T1 as has often been suggested. The differences in accuracy are shown to be significant because the evaluation at Tm using Eq. (35) is effectively a second-order accurate numerical method, with truncation error in the heat flux of O(ΔT3). In contrast, evaluation at T1 is a first-order approximation.
  3. (3)Equation (35) with εi(Tm) has excellent accuracy for |ΔT/T1| up to 20–30%. Beyond this range, the linearization error of Fig. 3 will cause increasing error for any surface.
  4. (4)Theory and examples for several nongray materials show that the gray-body approximation gives the wrong slope for heat flux as T2T1.
  5. (5)Calculations involving both the internal and external fractional functions can be conveniently implemented using the incomplete zeta function.
  6. (6)εi(Tm) should be used in calculating radiation thermal resistances for nongray surfaces, with Eq. (48).

Don Edwards was on the faculty of UCLA's CNTE Department when I studied there, although he was not my instructor. Nevertheless, I felt his positive influence through the curriculum and the other faculty in the heat transfer program.

  • c1 =

    first radiation constant, 3.741772 × 10−16 (W m2)

  • c2 =

    second radiation constant, 14387.77 (μm K)

  • eλ,b =

    blackbody monochromatic emissive power (W m−3)

  • f =

    (external) radiation fractional function

  • F =

    fi − f, Eq. (18)

  • fi =

    internal radiation fractional function

  • qnet =

    net radiant heat flux from surface 1 to 2 (W m−2)

  • re =

    electrical resistivity (Ω cm)

  • Rrad =

    radiation thermal resistance, Eq. (47) (K W−1)

  • T =

    temperature (K)

  • Tm =

    mean temperature, (T2 + T1)/2 (K)

  • T1 =

    temperature of surface (K)

  • T2 =

    temperature of environment (K)

  • X =

    c2/λT

  • Xz =

    finite value of X for which dF/dX = 0, 3.92069

  • Y =

    function in Appendix C, Eq. (D1)

  • αlwαsw =

    see Eq. (39)

  • α(T1, T2) =

    total hemispherical absorptivity

  • α(λ, T) =

    spectral hemispherical absorptivity

  • δT =

    T2 − T1 (K)

  • ΔT =

    T1 − T2 (K)

  • Δα =

    αlw − αsw

  • ε(T) =

    total hemispherical emissivity (external)

  • εn(T) =

    total normal emissivity (external)

  • εn(λT) =

    spectral normal emissivity

  • ε(λT) =

    spectral hemispherical emissivity

  • ε(θ,λ,T) =

    spectral directional emissivity

  • εi(T) =

    internal total hemispherical emissivity

  • εni(T) =

    internal total normal emissivity

  • ζ(s) =

    Riemann zeta function

  • ζ(Xs) =

    incomplete zeta function

  • θ =

    polar angle (rad)

  • λ =

    wavelength (m)

  • σ =

    Stefan–Boltzmann constant, 5.670367 × 10−8 (W m−2 K−4)

Appendix A: Incomplete Zeta Function

The radiation fractional function may be written in terms of the incomplete zeta function for convenience Display Formula

(A1)f(λT)=1σT40λ2πhco2λ5[exp(hco/kBTλ)1]dλ
Display Formula
(A2)=1σT42πkB4T4h3co2c2/λTt3et1dt

When λT, f =1 and the last equation yields the well-known result Display Formula

(A3)σT4=2πkB4T4h3co20t3et1dtζ(4)Γ(4)

where the Gamma function Γ(4) = 3! and the Riemann zeta function, ζ(4), has the indicated integral representation [14, Sec. 13.12]. A classical result due to Euler [15] gives ζ(4) = π4/90 (see also Ref. [16, Sec. 167]), from which we recover the usual definition of the Stefan–Boltzmann constant, σ. Returning to Eq. (A2) with this information, we have Display Formula

(A4)f(λT)=15π40t3et1dt15π40c2/λTt3et1dt
Display Formula
(A5)=115π40Xt3et1dt
Display Formula
(A6)=115π4Γ(4)ζ(X,4)

where X = c2/λT and we have identified the incomplete zeta function, ζ(X, s) [17, Sec. 8.22]. Hence Display Formula

(A7)f(λT)=190π4ζ(X,4)

Appendix B: Numerical Implementation

The integrals were computed using the GNU Scientific Library [18] with FFI bindings [19] to Lua code [20] under LuaLaTeX [21] using TEXShop over TEX Live [22]. The C code was compiled using XCode under Mac OS X. Integration was performed using GSL's QAG adaptive integration with 61 point Gauss–Kronrod rules. Convergence was checked, and the numerical values given in the text are believed to be accurate to the number of digits shown.

The incomplete zeta function was rendered in terms of the third-order Debye function, D3(X) = 3 ζ(X, 4) Γ(4)/X3, and computed using GSL's built-in routine. These routines were supplied to PGFPLOTS [23] to generate the charts.

The integrations of the discrete data for alumina were handled outside GSL, as described in the main text and Appendix E.

Readers interested in implementing these numerical methods may refer to the following tutorial publications and posts: [2427]. All of these packages (other than Mac OS X) are available online at no cost and are actively maintained.

Appendix C: The Constant Xz

Xz is the finite zero point of dF/dX, specifically the nonzero root of Display Formula

(C1)4(1eXz)=Xz

The equation may be solved by iteration to find Xz=3.92069. Alternatively, by defining y =4 − Xz, Eq. (C1) becomes Display Formula

(C2)yey=4e4

This equation may be solved using the Lambert W function Display Formula

(C3)Xz=4W(4e4)=3.92069

We may show by contradiction that Xz is irrational. Assuming Xz is rational, Xz = a/b for nonzero integers a and b. Then, the left-hand side of Eq. (C1) is 4 − 4ea/b; however, ea/b is itself irrational (e.g., see Ref. [28, Sec. 4.7]). Hence, 4 − 4ea/b is an irrational number and cannot equal the assumed right-hand side, a/b.

Diophantine approximations to Xz may be constructed by calculating a continued fraction representation of Xz through the usual process Display Formula

(C4)Xz=3.92069=3+11+111+11+11+1

Successive convergents (finite truncations) of the fraction give an infinite number of rational approximations to XzDisplay Formula

(C5)Xz{4,4712,,14938,24763,,1137290,}

The second of these is within 0.1% of the exact value; the last agrees to six digits.

Appendix D: Linearizations as Finite Difference Approximations

We may consider qnet to be a function of T that we wish to determine approximately at T = T2 by evaluating at other temperatures, i.e., T = T1 and T = Tm. With Display Formula

(D1)Y(T)=σε(T1)T14σα(T1,T)T4
we seek to estimate Display Formula
(D2)qnet=Y(T2)=σε(T1)T14σα(T1,T2)T24

A first-order Euler method starting at T = T1 would approximate Y (T2) based only on conditions at T1 as follows: Display Formula

(D3)Y(T2)Y(T1)+Y(T1)·δT=Y(T1)·δT
for δT = T2 − T1 (= −ΔT), which gives Eq. (7) or Eq. (14). The local truncation error of this single-step approximation is O(YΔT2) [29].3

A second-order Runge–Kutta method works from Tm with expansions toward both T1 and T2, subtracting the former from the latter Display Formula

(D4)Y(T2)=Y(Tm)+Y(Tm)δT2+Y(Tm)δT28+O(δT3)
Display Formula
(D5)Y(T1)=Y(Tm)Y(Tm)δT2+Y(Tm)δT28O(δT3)
Display Formula
(D6)Y(T2)=Y(T1)+Y(Tm)·δT+O(δT3)
Display Formula
(D7)Y(T2)Y(Tm)·δT

This single-step Runge–Kutta approximation has local truncation error O(YΔT3) [29]. In these calculations, α(λ, T1) is evaluated at the surface temperature, T1. With Eqs. (3) and (13)Display Formula

(D8)Y(Tm)=ddT(σT4α(T1,T))|Tm
Display Formula
(D9)=0α(λ,T1)eλ,bT|Tmdλ
Display Formula
(D10)=4σTm3·εi(Tm)
which results in Eq. (35). Given lower truncation error of the Runge–Kutta approach, we should expect significantly better accuracy using Eq. (35) than when using Eq. (7) or Eq. (14); and indeed the numerical results confirm this conclusion.

Appendix E: Integration of Directional Emissivity for Alumina

The spectral hemispherical emissivity was obtained by integration of the spectral directional emissivity data from Ref. [12] as Display Formula

(E1)ε(λ,T)=0π/2ε(θ,λ,T)sin(2θ)dθ

The spectral directional emissivity data are in 12 deg increments of polar angle θ from 0 deg to 72 deg. In all cases, the data are essentially constant from 0 to 36 deg, and this range was integrated analytically. From 36 deg to 84 deg, a five-point trapezoidal rule was used, and the integral from 84 deg to 90 deg was approximated as a trapezoid. The value at 90 deg was set to zero, in line with theory. This procedure was found to have a numerical truncation error of 1.0% for a gray surface.

The data for the lowest emissivities showed angular behavior consistent with a dielectric medium, as noted in Ref. [12], decreasing at higher angles. On this basis, the point at 84 deg was interpolated, using a value representative of variations at large angle for a dielectric: ε(84 deg, λ) ≈ 0.75 ε(72 deg, λ). Without large angle measurements or the complex refractive index, we cannot exclude the presence of the kind of peak emissivity above 80 deg predicted by Drude's model for metals [13]; and, in particular, the data at 17 and 23 μm do show a 30–50% increase in emissivity by 72 deg. Even so, a sensitivity analysis letting ε(84 deg, λ) ≈ 2.5 ε(72 deg, λ) increases the hemispherical emissivity by only about 5% of the previous estimate. Absent further data, there is not much basis upon which to adjust the calculation further.

Appendix F: Nondimensional Formulation of Model Surface Heat Exchange

Figures 4 and 5 plot heat exchange for a surface at T1 = Tz that switches between perfectly black and perfectly reflective behavior at wavelength λz, where Xz = c2/λzTz = 3.92069, for the case T1 = Tz = 300 K. The four equations used may be nondimensionalized as follows: Display Formula

(F1)exact:qnetσT14=ε(Tz)α(Tz,T2)(T2T1)4
Display Formula
(F2)gray:qnetσT14ε(Tz)[1(T2T1)4]
Display Formula
(F3)εi(Tz):qnetσT144εi(Tz)ΔTT1
Display Formula
(F4)εi(Tm):qnetσT144εi(Tm)(TmT1)3ΔTT1

In all cases, qnet/σT14=fn(T2/T1), since Tm/T1 and ΔT/T1 may both be written in terms of T2/T1.

Surface Black Below λz.

From Eqs. (27), (45), and (46) with αsw = 1, Δα = −1, and λc = λz, we have Display Formula

(F5)ε(Tz)=0.4177andεi(Tz)=0.6017
Display Formula
(F6)εi(Tm)=1[90π4ζ(Xc,m,4)F(Xc,m)]
Display Formula
(F7)α(Tz,T2)=190π4ζ(Xc,2,4)

where Xc,m=c2/λzTm=Xz(T1/Tm)=3.92069(T1/Tm) and Xc,2=c2/λzT2=Xz(T1/T2)=3.92069(T1/T2).

The curves for this case are plotted in Fig. 10. Accuracy of the linearization-based εi(Tm) is excellent over the range T2/T1 ≈ 1 ± 0.30 in even this extreme case of nongray behavior. The gray approximation strongly overpredicts qnet when T2 > T1 because it underestimates the absorption of irradiation.

Surface Black Above λz.

From Eqs. (28), (45), and (46) with αsw = 0 and Δα = 1, we have Display Formula

(F8)ε(Tz)=0.5823andεi(Tz)=0.3983
Display Formula
(F9)εi(Tm)=[90π4ζ(Xc,m,4)F(Xc,m)]
Display Formula
(F10)α(Tz,T2)=90π4ζ(Xc,2,4)

The curves for this case are plotted in Fig. 11. As before, the linearization-based εi(Tm) is very accurate over the range T2/T1 = 1 ± 0.30. The gray approximation strongly underpredicts qnet when T2 > T1 because it overestimates the absorption of irradiation.

The internal emissivity at Tm is shown for both surfaces in Fig. 12.

Edwards, D. K. , 1969, “ Radiative Transfer Characteristics of Materials,” ASME J. Heat Transfer, 91(1), pp. 1–15. [CrossRef]
Czerny, M. , and Walther, A. , 1961, Tables of the Fractional Functions of the Planck Radiation Law, Springer-Verlag, Berlin.
Edwards, D. K. , and Nelson, K. E. , 1961, “ Maximum Error in Total Emissivity Measurements Due to Non-Grayness of Samples,” ARS J., 31(7), pp. 1021–1022.
Edwards, D. K. , 1970, “ Thermal Radiation Measurements,” Measurement Techniques in Heat Transfer (AGARDograph, Vol. 130), E. R. G. Eckert and R. J. Goldstein , eds., Advisory Group for Aerospace Research and Development of NATO, Technivision Services, Slough, UK, Chap. 9.
Edwards, D. K. , 1981, Radiation Heat Transfer Notes, Hemisphere Publishing, New York.
Edwards, D. K. , Denny, V. E. , and Mills, A. F. , 1979, Transfer Processes, 2nd ed., Hemisphere Publishing, Washington, DC.
Mills, A. F. , 1999, Heat Transfer, 2nd ed., Prentice Hall, Upper Saddle River, NJ.
Rosseland, S. , 1936, Theoretical Astrophysics, Oxford University Press, London.
Brewster, M. Q. , 1992, Thermal Radiative Transfer and Properties, Wiley, New York.
Modest, M. F. , 1993, Radiative Heat Transfer, McGraw-Hill, New York.
Lienhard , J. H., IV , and Lienhard, J. H., V , 2018, A Heat Transfer Textbook, 4th ed., Phlogiston Press, Cambridge, MA.
Teodorescu, G. , and Jones, P. D. , 2008, “ Spectral and Directional Emittance of Alumina at 823 K,” J. Mater. Sci., 43(22), pp. 7225–7229. [CrossRef]
Baehr, H. D. , and Stephan, K. , 1998, Heat and Mass Transfer, Springer, Berlin.
Whittaker, E. T. , and Watson, E. N. , 1920, A Course of Modern Analysis, 3rd ed., Cambridge University Press, Cambridge, UK.
Euler, L. , 1740, “ De summis serierum reciprocarum,” Commentarii Academiae Scientiarum Petropolitanae, Vol. 7, pp. 123–134, First communicated to Daniel Bernoulli in 1734 and read before the St. Petersburg Academy in December 1735.
Euler, L. , 1748, Introductio in analysin infinitorum, Vol. 1, Apud Marcum-Michaelem Bousquet & Socios, Lausannae, Switzerland.
Olver, F. W. J. , Olde Daalhuis, A. B. , Lozier, D. W. , Schneider, B. I. , Boisvert, R. F. , Clark, C. W. , Miller, B. R. , and Saunders, B. V. , eds., 2017, “NIST Digital Library of Mathematical Functions, Release 1.0.17,” National Institute of Standards and Technology, Gaithersburg, MD, accessed Dec. 22, 2018 , http://dlmf.nist.gov/
Galassi, M. , Davies, J. , Theiler, J. , Gough, B. , Jungman, G. , Alken, P. , Booth, M. , Rossi, F. , and Ulerich, R. , 2017, “ GNU Scientific Library Release 2.4,” Free Software Foundation, Boston, MA, accessed Dec. 12, 2018, https://www.gnu.org/software/gsl/
Pall, M., 2017, “FFI Library,” luajit, Muenchen, Germany, accessed Mar. 1, 2018, http://www.luatex.org/svn/branches/0.79/source/libs/luajit/LuaJIT-2.0.3/doc/ext_ffi.html
Pontifical Catholic University of Rio de Janeiro, 2017, “ Lua 5.3 Reference Manual,” Pontifical Catholic University of Rio de Janeiro, Rio de Janeiro, Brazil, accessed June 18, 2018, https://www.lua.org/manual/5.3/
Pégourié-Gonnard, M. , 2013, “ A Guide to LuaLaTeX,” DANTE eV, Heidelberg, Germany, accessed May 5, 2013, http://dante.ctan.org/tex-archive/info/luatex/lualatex-doc/lualatex-doc.pdf
Koch, R. , 2018, “ TEXShop, Version 4.01,” University of Oregon, Eugene, OR, accessed Apr. 16, 2018, http://pages.uoregon.edu/koch/texshop/obtaining.html
Feuersänger, C. , 2017, “ Manual for Package PGFPLOTS, Version 1.15,” SourceForge Media, La Jolla, CA, accessed June 5, 2017, http://sourceforge.net/projects/pgfplots
Montijano, J. I. , Pérez, M. , Rández, L. , and Varona, J. L. , 2014, “ Numerical Methods With LuaLaTeX,” TUGboat, 35(1), pp. 51–56. https://tug.org/TUGboat/tb35-1/tb109montijano.pdf
Menke, H. , 2017, “ Plot Complete Fermi-Dirac Integral in LuaLaTeX,” Stack Overflow, New York, Nov. 29, 2017, https://tex.stackexchange.com/questions/403766/plot-complete-fermi-dirac-integral-in-lualatex
Menke, H. , 2018, “ Tutorial: Using External C Libraries With the LuaLaTeX FFI,” TUGboat, 39(1), pp. 37–40. https://www.tug.org/TUGboat/tb39-1/tb121menke-ffi.pdf
Dyer, R. , 2015, “ Compiling the GSL Library for OSX,” Virginia Commonwealth University, Richmond, VA, accessed July 29, 2015, https://dyerlab.ces.vcu.edu/2015/07/29/compiling-the-gsl-library-for-osx/
Hardy, G. H. , and Wright, E. M. , 1979, An Introduction to the Theory of Numbers, 5th ed., Oxford University Press, Oxford, UK.
Press, W. H. , Teukolsky, S. A. , Vetterling, W. T. , and Flannery, B. P. , 1988, Numerical Recipes in C, Cambridge University Press, Cambridge, UK.
Hornbeck, R. W. , 1975, Numerical Methods, Quantum Publishers, New York.
Copyright © 2019 by ASME; reuse license CC-BY 4.0
View article in PDF format.

References

Edwards, D. K. , 1969, “ Radiative Transfer Characteristics of Materials,” ASME J. Heat Transfer, 91(1), pp. 1–15. [CrossRef]
Czerny, M. , and Walther, A. , 1961, Tables of the Fractional Functions of the Planck Radiation Law, Springer-Verlag, Berlin.
Edwards, D. K. , and Nelson, K. E. , 1961, “ Maximum Error in Total Emissivity Measurements Due to Non-Grayness of Samples,” ARS J., 31(7), pp. 1021–1022.
Edwards, D. K. , 1970, “ Thermal Radiation Measurements,” Measurement Techniques in Heat Transfer (AGARDograph, Vol. 130), E. R. G. Eckert and R. J. Goldstein , eds., Advisory Group for Aerospace Research and Development of NATO, Technivision Services, Slough, UK, Chap. 9.
Edwards, D. K. , 1981, Radiation Heat Transfer Notes, Hemisphere Publishing, New York.
Edwards, D. K. , Denny, V. E. , and Mills, A. F. , 1979, Transfer Processes, 2nd ed., Hemisphere Publishing, Washington, DC.
Mills, A. F. , 1999, Heat Transfer, 2nd ed., Prentice Hall, Upper Saddle River, NJ.
Rosseland, S. , 1936, Theoretical Astrophysics, Oxford University Press, London.
Brewster, M. Q. , 1992, Thermal Radiative Transfer and Properties, Wiley, New York.
Modest, M. F. , 1993, Radiative Heat Transfer, McGraw-Hill, New York.
Lienhard , J. H., IV , and Lienhard, J. H., V , 2018, A Heat Transfer Textbook, 4th ed., Phlogiston Press, Cambridge, MA.
Teodorescu, G. , and Jones, P. D. , 2008, “ Spectral and Directional Emittance of Alumina at 823 K,” J. Mater. Sci., 43(22), pp. 7225–7229. [CrossRef]
Baehr, H. D. , and Stephan, K. , 1998, Heat and Mass Transfer, Springer, Berlin.
Whittaker, E. T. , and Watson, E. N. , 1920, A Course of Modern Analysis, 3rd ed., Cambridge University Press, Cambridge, UK.
Euler, L. , 1740, “ De summis serierum reciprocarum,” Commentarii Academiae Scientiarum Petropolitanae, Vol. 7, pp. 123–134, First communicated to Daniel Bernoulli in 1734 and read before the St. Petersburg Academy in December 1735.
Euler, L. , 1748, Introductio in analysin infinitorum, Vol. 1, Apud Marcum-Michaelem Bousquet & Socios, Lausannae, Switzerland.
Olver, F. W. J. , Olde Daalhuis, A. B. , Lozier, D. W. , Schneider, B. I. , Boisvert, R. F. , Clark, C. W. , Miller, B. R. , and Saunders, B. V. , eds., 2017, “NIST Digital Library of Mathematical Functions, Release 1.0.17,” National Institute of Standards and Technology, Gaithersburg, MD, accessed Dec. 22, 2018 , http://dlmf.nist.gov/
Galassi, M. , Davies, J. , Theiler, J. , Gough, B. , Jungman, G. , Alken, P. , Booth, M. , Rossi, F. , and Ulerich, R. , 2017, “ GNU Scientific Library Release 2.4,” Free Software Foundation, Boston, MA, accessed Dec. 12, 2018, https://www.gnu.org/software/gsl/
Pall, M., 2017, “FFI Library,” luajit, Muenchen, Germany, accessed Mar. 1, 2018, http://www.luatex.org/svn/branches/0.79/source/libs/luajit/LuaJIT-2.0.3/doc/ext_ffi.html
Pontifical Catholic University of Rio de Janeiro, 2017, “ Lua 5.3 Reference Manual,” Pontifical Catholic University of Rio de Janeiro, Rio de Janeiro, Brazil, accessed June 18, 2018, https://www.lua.org/manual/5.3/
Pégourié-Gonnard, M. , 2013, “ A Guide to LuaLaTeX,” DANTE eV, Heidelberg, Germany, accessed May 5, 2013, http://dante.ctan.org/tex-archive/info/luatex/lualatex-doc/lualatex-doc.pdf
Koch, R. , 2018, “ TEXShop, Version 4.01,” University of Oregon, Eugene, OR, accessed Apr. 16, 2018, http://pages.uoregon.edu/koch/texshop/obtaining.html
Feuersänger, C. , 2017, “ Manual for Package PGFPLOTS, Version 1.15,” SourceForge Media, La Jolla, CA, accessed June 5, 2017, http://sourceforge.net/projects/pgfplots
Montijano, J. I. , Pérez, M. , Rández, L. , and Varona, J. L. , 2014, “ Numerical Methods With LuaLaTeX,” TUGboat, 35(1), pp. 51–56. https://tug.org/TUGboat/tb35-1/tb109montijano.pdf
Menke, H. , 2017, “ Plot Complete Fermi-Dirac Integral in LuaLaTeX,” Stack Overflow, New York, Nov. 29, 2017, https://tex.stackexchange.com/questions/403766/plot-complete-fermi-dirac-integral-in-lualatex
Menke, H. , 2018, “ Tutorial: Using External C Libraries With the LuaLaTeX FFI,” TUGboat, 39(1), pp. 37–40. https://www.tug.org/TUGboat/tb39-1/tb121menke-ffi.pdf
Dyer, R. , 2015, “ Compiling the GSL Library for OSX,” Virginia Commonwealth University, Richmond, VA, accessed July 29, 2015, https://dyerlab.ces.vcu.edu/2015/07/29/compiling-the-gsl-library-for-osx/
Hardy, G. H. , and Wright, E. M. , 1979, An Introduction to the Theory of Numbers, 5th ed., Oxford University Press, Oxford, UK.
Press, W. H. , Teukolsky, S. A. , Vetterling, W. T. , and Flannery, B. P. , 1988, Numerical Recipes in C, Cambridge University Press, Cambridge, UK.
Hornbeck, R. W. , 1975, Numerical Methods, Quantum Publishers, New York.

Figures

Grahic Jump Location
Fig. 1

Internal and external radiation fractional functions and their difference

Grahic Jump Location
Fig. 2

Total absorptivity for both of the example surfaces, α(300 K, T2)

Grahic Jump Location
Fig. 3

Linearization error for a black surface about T1 or aboutTm

Grahic Jump Location
Fig. 4

Comparison of models for qnet (300 K surface, black below 12.23 μm)

Grahic Jump Location
Fig. 5

Comparison of models for qnet (300 K surface, black above 12.23 μm)

Grahic Jump Location
Fig. 6

Spectral normal emissivity of alumina at 823 K (Reprinted by permission from Teodorescu and Jones [12]. Copyright 2008 by Springer.)

Grahic Jump Location
Fig. 7

Comparison of models for alumina surface at T1 = 823 K

Grahic Jump Location
Fig. 8

Comparison of models for platinum surface at T1 = 373 K

Grahic Jump Location
Fig. 9

Model absorptivity for solar reflector, based on soft anodized aluminum

Grahic Jump Location
Fig. 10

Models for qnet (black below λz)

Grahic Jump Location
Fig. 11

Models for qnet (black above λz)

Grahic Jump Location
Fig. 12

Internal emissivity, εi(Tm), for model surfaces

Tables

Table Grahic Jump Location
Table 1 Selective solar reflector: soft anodized aluminum at T1 = 360 K with T2 = 290 K. Heat flux in W m−2.

Errata

Some tools below are only available to our subscribers or users with an online account.

Related Content

Customize your page view by dragging and repositioning the boxes below.

Related Journal Articles
Related eBook Content
Topic Collections

Sorry! You do not have access to this content. For assistance or to subscribe, please contact us:

  • TELEPHONE: 1-800-843-2763 (Toll-free in the USA)
  • EMAIL: asmedigitalcollection@asme.org
Sign In