Analytical Solution to Transient Asymmetric Heat ... - Semantic Scholar
Prashant K. Jain1 Suneet Singh Rizwan-uddin Department of Nuclear, Plasma and Radiological Engineering, University of Illinois at Urbana–Champaign, 216 Talbot Lab, 104 South Wright Street, Urbana, IL 61801
Analytical Solution to Transient Asymmetric Heat Conduction in a Multilayer Annulus In this paper, we present an analytical double-series solution for the time-dependent asymmetric heat conduction in a multilayer annulus. In general, analytical solutions in multidimensional Cartesian or cylindrical 共r , z兲 coordinates suffer from existence of imaginary eigenvalues and thus may lead to numerical difficulties in computing analytical solution. In contrast, the proposed analytical solution in polar coordinates (2D cylindrical) is “free” from such imaginary eigenvalues. Real eigenvalues are obtained by virtue of precluded explicit dependence of transverse (radial) eigenvalues on those in the other direction. 关DOI: 10.1115/1.2977553兴 Keywords: heat conduction, layered annulus, analytical method
1
Introduction
In modern engineering applications, multilayer components are extensively used due to the added advantage of combining physical, mechanical, and thermal properties of different materials. Many of these applications require a detailed knowledge of transient temperature and heat-flux distribution within the component layers. Both analytical and numerical techniques may be used to solve such problems. Nonetheless, numerical solutions are preferred and prevalent in practice, due to either unavailability or higher mathematical complexity of the corresponding exact solutions. Rather limited use of analytical solutions should not diminish their merit over numerical ones; since exact solutions, if available, provide an insight into the governing physics of the problem, which is typically missing in any numerical solution. Moreover, analyzing closed-form solutions to obtain optimal design options for any particular application of interest is relatively simpler. In addition, exact solutions find their applications in validating and comparing various numerical algorithms to help improve computational efficiency of computer codes that currently rely on numerical techniques. Although multilayer heat conduction problems have been studied in great detail and various solution methods—including orthogonal and quasiorthogonal expansion technique, Laplace transform method, Green’s function approach, finite integral transform technique 关1–11兴—are readily available; there is a continued need to develop and explore novel methods to solve problems for which exact solutions still do not exist. One such problem is to determine exact unsteady temperature distribution in polar coordinates 共r , 兲 with multiple layers in the radial direction. Salt 关12,13兴 addressed time-dependent heat conduction problem by orthogonal expansion technique, in a two-dimensional composite slab 共Cartesian geometry兲 with no internal heat source, subjected to homogenous boundary conditions. Later, Mikhailov and Ozisik 关14兴 solved the 3D transient conduction problem in a Cartesian nonhomogenous finite medium. More recently, Haji-Sheikh and Beck 关15兴 applied Green’s function approach to develop transient temperature solutions in a 3D Cartesian two-layer orthotropic medium including the effects of contact resistance. Lu et al. 关5兴 developed a novel method by combining Laplace transform 1 Corresponding author. Contributed by the Heat Transfer Division of ASME for publication in the JOURNAL OF HEAT TRANSFER. Manuscript received February 12, 2008; final manuscript received April 21, 2008; published online October 20, 2008. Review conducted by Peter Vadasz.
Journal of Heat Transfer
method and separation of variables method to solve multidimensional transient heat conduction problem in a rectangular and cylindrical multilayer slab with time-dependent periodic boundary condition. The treatment in the cylindrical coordinates is, however, restricted to the r – z coordinates. Eigenfunction expansion method is applied by de Monte 关16兴 to solve the unsteady heat conduction problem in a two-dimensional two-layer isotropic slab subjected to homogenous boundary conditions. Feng-Bin Yeh 关17兴 applied the method of separation of variables to solve plasma heating of a one-dimensional two-layer composite slab with layers in imperfect thermal contact. The brief review of relevant literature is by no means exhaustive. However, a literature survey showed that the analytical solution for unsteady temperature distribution in a multilayer annular geometry has not yet been developed. Recently, an exact solution based on the separation of variables method is developed by Singh et al. 关18兴 for multilayer heat conduction in polar coordinates. However, that exact solution is applicable only to domains with pie slice geometry 共 ⬍ 2, where is the angle subtended by the layers兲. Numerous applications involving multilayer cylindrical geometry require evaluation of temperature distribution in complete disk-type 共i.e., = 2兲 layers. One typical example is a nuclear fuel rod, which consists of concentric layers of different materials and often subjected to asymmetric boundary conditions. Moreover, several other applications including multilayer insulation materials, double heat-flux conductimeter, typical laser absorption calorimetry experiments, cryogenic systems, and other cylindrical building structures would benefit from such analytical solutions. This paper extends the solution approach developed by Singh et al. 关18兴 for such applications and presents an analytical double-series solution for the time-dependent asymmetric heat conduction in a multilayer annulus. Solution is valid for any combinations of time-independent, inhomogeneous boundary conditions at the inner and outer radii of the domain. The results for an illustrative problem involving a three-layer annulus subjected to asymmetric heat-flux are also presented.
2
Mathematical Formulation
Consider an n-layer annulus 共r0 艋 r 艋 rn兲, as shown schematically in Fig. 1. All the layers are assumed to be isotropic in thermal properties and are in perfect thermal contact. Let ki and ␣i be the temperature independent thermal conductivity and thermal diffusivity of the ith layer. At t = 0, each ith layer is at a specified temperature f i共r , 兲 and time-independent heat sources gi共r , 兲 are switched on for t ⬎ 0. Both the inner 共i = 1, r = r0兲 as well as the
Copyright © 2009 by ASME
JANUARY 2009, Vol. 131 / 011304-1
Downloaded 01 Mar 2010 to 130.126.32.13. Redistribution subject to ASME license or copyright; see http://www.asme.org/terms/Terms_Use.cfm
ki
Ti Ti−1 共ri−1, ,t兲 = ki−1 共ri−1, ,t兲 r r
共7兲
Initial condition 共for ri−1 艋 r 艋 ri and 0 艋 艋 2, where i = 1 , 2 , . . . , n兲: Ti共r, ,t = 0兲 = f i共r, 兲
共8兲
Boundary conditions either of the first, second, or third kind may be imposed at r = r0 and r = rn by choosing the appropriate coefficients in Eqs. 共2兲 and 共3兲. However, the case in which Bin and Bout are simultaneously zero is not considered. In addition, asymmetric boundary conditions can be applied by choosing -dependent Cin and Cout. Furthermore, multiple layers with zero inner radius 共r0 = 0兲 can be simulated by assigning zero values to constants Bin and Cin in Eq. 共2兲. It should be noted that the formulation presented in this paper only applies toward timeindependent boundary conditions and/or source terms due to the limitation of the separation of variables method. This solution methodology cannot be extended to include the effects of timedependence in boundary conditions and/or sources. Such problems can be solved analytically using the finite integral transform technique 关20,21兴. Fig. 1 Schematic of an n-layer annulus. ith layer has an inner and outer radii equal to ri−1 and ri, respectively.
outer 共i = n, r = rn兲 surfaces of the annulus may be subjected to any combination of temperature and heat-flux boundary conditions. 共Since perfect thermal contact between the adjacent layers is seldom observed in real materials, dealing with imperfect contact would require explicit modeling of the thermal resistance at the layer interfaces 关11,17,19兴. For such cases, the temperature at the contact interfaces will not be continuous.兲 Under these assumptions, the governing heat conduction equation along with the boundary and initial conditions are as follows. Governing equation 共for ri−1 艋 r 艋 ri, 0 艋 艋 2, and t ⬎ 0, where i = 1 , 2 , . . . , n兲:
冉
冊
3
Solution Methodology
In order to apply the separation of variables method, which is only applicable to homogenous problems, the nonhomogenous problem has to be split 关21兴 into: 共1兲 homogenous transient problem, and 共2兲 nonhomogenous steady-state problem. This is accomplished by splitting Ti共r , , t兲 in the governing Eqs. 共1兲–共8兲 as ¯T 共r , , t兲 + T 共r , 兲, where ¯T 共r , , t兲 is the “complementary” i
ss,i
3.1 Homogenous Transient Problem. Homogenized complementary transient equations corresponding to Eqs. 共1兲–共8兲 are as follows. Governing equation 共for ri−1 艋 r 艋 ri, 0 艋 艋 2 and t ⬎ 0, where i = 1 , 2 , . . . , n兲:
冉
1 Ti 1 2T i gi共r, 兲 1 Ti 共r, ,t兲 = r 共r, ,t兲 + 2 2 共r, ,t兲 + ␣i t r r r r ki 共1兲 Boundary conditions: •
•
Tn 共rn, ,t兲 + BoutTn共rn, ,t兲 = Cout共兲 Aout r •
•
Inner surface of the first layer 共for 0 艋 艋 2 and t ⬎ 0兲 Ain
共2兲 •
•
¯T1 共r0, ,t兲 + Bin¯T1共r0, ,t兲 = 0 r
共10兲
Outer surface of the nth layer 共for 0 艋 艋 2 and t ⬎ 0兲 Aout
共3兲
Periodic boundary conditions 共for ri−1 艋 r 艋 ri and t ⬎ 0, where i = 1 , 2 , . . . , n兲
共9兲
Boundary conditions: •
Outer surface of the nth layer 共for 0 艋 艋 2 and t ⬎ 0兲
冊
1 ¯Ti 1 2¯Ti 1 ¯Ti 共r, ,t兲 = r 共r, ,t兲 + 2 2 共r, ,t兲 ␣i t r r r r
Inner surface of the first layer 共for 0 艋 艋 2 and t ⬎ 0兲
T1 共r0, ,t兲 + BinT1共r0, ,t兲 = Cin共兲 Ain r
i
transient part and Tss,i共r , 兲 is the steady-state part of the solution.
¯Tn 共rn, ,t兲 + Bout¯Tn共rn, ,t兲 = 0 r
共11兲
Periodic boundary conditions 共for ri−1 艋 r 艋 ri and t ⬎ 0, where i = 1 , 2 , . . . , n兲
Ti共r, = 0,t兲 = Ti共r, = 2,t兲
共4兲
¯T 共r, = 0,t兲 = ¯T 共r, = 2,t兲 i i
共12兲
Ti Ti 共r, = 0,t兲 = 共r, = 2,t兲
共5兲
¯Ti ¯Ti 共r, = 0,t兲 = 共r, = 2,t兲
共13兲
Interface of the 共i − 1兲st and the ith layer 共for 0 艋 艋 2 and t ⬎ 0, where i = 2 , . . . , n兲 Ti共ri−1, ,t兲 = Ti−1共ri−1, ,t兲
011304-2 / Vol. 131, JANUARY 2009
共6兲
•
Interface of the 共i − 1兲st and the ith layer 共for 0 艋 艋 2 and t ⬎ 0 where i = 2 , . . . , n兲 ¯T 共r , ,t兲 = ¯T 共r , ,t兲 i i−1 i−1 i−1
共14兲
Transactions of the ASME
Downloaded 01 Mar 2010 to 130.126.32.13. Redistribution subject to ASME license or copyright; see http://www.asme.org/terms/Terms_Use.cfm
ki
¯Ti ¯Ti−1 共ri−1, ,t兲 = ki−1 共ri−1, ,t兲 r r
⬁
¯T 共r, ,t兲 = i
共15兲
冊
1 Tss,i 1 2Tss,i gi共r, 兲 共r, 兲 + 2 r 共r, 兲 + = 0, r r r r 2 ki ri−1 艋 r 艋 ri, 1 艋 i 艋 n
+
共17兲
Tss,1 共r0, 兲 + BinTss,1共r0, 兲 = Cin共兲 r
Tss,n Aout 共rn, 兲 + BoutTss,n共rn, 兲 = Cout共兲 r
共19兲
Tss,i Tss,i 共r, = 0兲 = 共r, = 2兲
共21兲
Tss,i Tss,i−1 共ri−1, 兲 = ki−1 共ri−1, 兲 r r
共22兲 共23兲
Solution to the Homogenous Transient Problem 4.1
¯T 共r, ,t兲 = R 共r兲⌰共兲⌫ 共t兲 i i i
冊冉
冊
1 d dRi共r兲 m r + i2 − 2 Ri共r兲 = 0 dr r dr r 2
in ri−1 艋 r 艋 ri, 共25兲
where i = 1,2, . . . ,n d2⌰共兲 + m2⌰共兲 = 0 d2 1 d⌫i共t兲 + i2⌫i共t兲 = 0 ␣i dt
in 0 艋 艋 2
Rimp共impr兲sin共m兲
共28兲
共29兲
i = 1,2, . . . ,n
The radial 共transverse兲 eigenfunction, Rimp共impr兲 in Eq. 共28兲 is
for t ⬎ 0,
兺
ki ␣i
共30兲
冕
ri
rRimp共impr兲Rimq共imqr兲dr =
ri−1
冋
0
if p ⫽ q
Nrmp if p = q
册
where Jm and Y m are Bessel functions of the first and second kind of order m, respectively. Nrmp is normalization integral in the r-direction. For the angular eigenfunctions ⌰m共兲—formed via combination of constant sin共m兲 and cos共m兲—the standard orthogonality condition is valid 关21兴. 4.3 Radial (Transverse) Eigencondition. Application of the boundary conditions 共Eqs. 共10兲 and 共11兲兲 and interface conditions 共Eqs. 共14兲 and 共15兲兲 to the transverse eigenfunction Rimp共impr兲 yields a 共2n ⫻ 2n兲 matrix for each integer value of m. Transverse eigencondition is obtained by setting the determinant of this matrix equal to zero. Roots of which, in turn, yield the infinite number of eigenvalues 1mp corresponding to the first layer for each integer value of m. 4.4 Determination of Coefficients aimp and bimp. Coefficients aimp and bimp in the radial eigenfunction Rimp共impr兲 共see Eq. 共30兲兲 are evaluated from the following recurrence relationship, obtained from the ith interface condition, valid for i 苸 共1 , n − 1兲:
冉 冊冉 bi+1,mp
=
Jm共i+1,mpri兲
冊 冊冉 冊
Y m共i+1,mpri兲
−1
⬘ 共i+1,mpri兲 ki+1Y m⬘ 共i+1,mpri兲 ki+1Jm
⫻
冉
Jm共impri兲
Y m共impri兲
⬘ 共impri兲 kiY m⬘ 共impri兲 k iJ m
aimp
bimp
共32兲
where b1mp = −共C1in / C2in兲a1mp, and a1mp is arbitrary. 4.5 Determination of Coefficients D0p, Dmp, and Emp. Coefficients D0p, Dmp, and Emp in Eq. 共28兲 are evaluated by applying the initial condition 共Eq. 共16兲兲 and then making use of the orthogonality conditions in the radial and angular directions, as follows n
D0p =
兺
1 ki 2Nr0p i=1 ␣i
冕冕 2
0
ri
¯ 共r, ,t = 0兲drd rRi0p共i0pr兲T i
ri−1
where i = 1,2, . . . ,n
共33兲
are constants of separation.
4.2 General Solution. In view of the ODEs listed above, a general solution for Eq. 共9兲 may be realized as Journal of Heat Transfer
n
共26兲
共27兲 where
2 −␣iimp t
where continuity of the heat-flux at the layer interfaces requires the following relationship between the ith imp and 1mp to hold 关15,16,21,22兴,
共24兲
in Eq. 共9兲 and applying separation of variables yield the following ordinary differential equations 共ODEs兲.
2i
mpe
m=1 p=1
ai+1,mp
Separation of Variables. Substituting the product form
冉
兺兺E
共31兲
共20兲
ki
Rimp共impr兲cos共m兲
⬁
⬁
+
i=1
Tss,i共r, = 0兲 = Tss,i共r, = 2兲
Tss,i共ri−1, 兲 = Tss,i−1共ri−1, 兲
4
共18兲
Interface of the 共i − 1兲st and the ith layer 共for 0 艋 艋 2, where i = 2 , . . . , n兲
•
2 −␣iimp t
and the corresponding orthogonality condition is 关18兴
Periodic boundary conditions 共for ri−1 艋 r 艋 ri, where i = 1 , 2 , . . . , n兲
•
mpe
Rimp共impr兲 = aimpJm共impr兲 + bimpY m共impr兲
Outer surface of the nth layer 共for 0 艋 艋 2兲
•
兺兺D
imp = 1mp冑␣1/␣i
Inner surface of the first layer 共for 0 艋 艋 2兲 Ain
Ri0p共i0pr兲
m=1 p=1
Boundary conditions: •
2 −␣ii0p t
⬁
⬁
共16兲
3.2 Inhomogeneous Steady-State Problem. Inhomogeneous steady-state equations corresponding to Eqs. 共1兲–共8兲 are as follows. Governing equation 共for ri−1 艋 r 艋 ri and 0 艋 艋 2, where i = 1 , 2 , . . . , n兲:
冉
0pe
p=1
Initial condition 共for ri−1 艋 r 艋 ri and 0 艋 艋 2, where i = 1 , 2 , . . . , n兲: ¯T 共r, ,t = 0兲 = f 共r, 兲 − T 共r, 兲 i i ss,i
兺D
n
兺
1 ki Dmp = Nrmp i=1 ␣i
冕冕 2
0
¯ 共r, ,t = 0兲drd ⫻T i
ri
rRimp共impr兲cos共m兲
ri−1
共34兲
JANUARY 2009, Vol. 131 / 011304-3
Downloaded 01 Mar 2010 to 130.126.32.13. Redistribution subject to ASME license or copyright; see http://www.asme.org/terms/Terms_Use.cfm
n
Emp =
兺
1 ki Nrmp i=1 ␣i
冕冕 2
0
ri
¯ 共r, ,t = 0兲drd rRimp共impr兲sin共m兲T i
ri−1
共35兲
5
Absence of imaginary Radial Eigenvalues
In general, transverse eigenvalues for multilayer timedependent heat conduction problems in Cartesian coordinates may be imaginary 关16兴. Same is true for 2D 共r , z兲 cylindrical coordinates. The eigenvalues are imaginary due to the explicit dependence of the transverse eigenvalues on those in the remaining direction共s兲. For example, in a two-layer homogenous heat conduction problem with layers in the x-direction, the general solution in the ith layer is 关15,16兴 ⬁
⬁
Ti共x,y,t兲 =
兺兺Z
mpe
2 2 −␣i共imp +m 兲t
Ximp共impx兲Y m共my兲 共36兲
m=1 p=1
where imp and m are eigenvalues in the x- and y-directions, respectively; and Ximp 共impx兲 and Y m共my兲 are the corresponding eigenfunctions. For heat-flux to be continuous at the interface ∀t 2 2 2 2 ␣1共1mp + m 兲 = ␣2共2mp + m 兲
which implies
2mp =
冑
共37兲
冉 冊
␣1 2 ␣1 2 + − 1 m ␣2 1mp ␣2
共38兲
Clearly, the relationship above may result in either real or imaginary transverse eigenvalues 关13,15,16兴. However, in the present case, similar considerations led to Eq. 共29兲, which is similar to what has been established for 1D multilayer time-dependent problems and eliminates the possibility of imaginary eigenvalues. Moreover, physical considerations dictate that eigenvalues should be real 共both in the present and the corresponding 1D case兲 because imaginary eigenvalues will result in exponentially growing temperatures. In contrast, the solution given in Eq. 共36兲 can have a physically realizable solution even in the case where eigenvalues in one of the directions are imaginary. It should be noted that, though there is no explicit dependence between radial and angular eigenvalues, the order of the Bessel functions constituting radial eigenfunctions is determined by the angular eigenvalues. Hence, the radial eigenvalues implicitly depend on the angular eigenvalues. Moreover, unlike in Cartesian coordinates, this implicit dependence does not vanish even if ␣1 = ␣i共i ⫽ 1兲. In fact, it exists even for single-layer problems.
Fig. 2 Asymmetric heat conduction in a three-layer annulus. Each layer has a different thermal conductivity „ki… and thermal diffusivity „␣i…. The lower-half of the annulus „ Ï Ï 2… is kept insulated, while the upper-half „0 Ï Ï … is subjected to a -dependent incoming heat-flux.
冉
冊
1 d dTˆs,im共r兲 m2 gˆs,im共r兲 r − 2 Tˆs,im共r兲 + =0 dr r dr r ki
冉
冊
1 d dTˆi0共r兲 gˆi0共r兲 r + =0 dr r dr ki
共41兲
共42兲
where the source term gi共r , 兲 is expanded in a generalized Fourier series as ⬁
gi共r, 兲 = gˆi0共r兲 +
兺
⬁
gˆc,im共r兲cos共m兲 +
m=1
兺 gˆ
s,im共r兲sin共m兲
m=1
共43兲 and
6 Solution to the Inhomogeneous Steady-State Problem The inhomogeneous steady-state problem is solved using the eigenfunction expansion method. The steady-state temperature distribution, governed by Eq. 共17兲, may be written as a generalized Fourier series in terms of angular eigenfunctions, ⬁
Tss,i共r, 兲 = Tˆi0共r兲 +
兺 Tˆ
s,im共r兲sin共m兲,
ri−1 艋 r 艋 ri,
1艋i艋n
m=1
共39兲 Substituting Eq. 共39兲 in Eq. 共17兲 leads to the following ODEs:
冉
冕
1
冕
gˆi0共r兲 ⬅
⬁
兺 Tˆ
gˆs,im共r兲 ⬅
1
c,im共r兲cos共m兲
m=1
+
gˆc,im共r兲 ⬅
冊
1 d dTˆc,im共r兲 m2 gˆc,im共r兲 r − 2 Tˆc,im共r兲 + =0 dr r dr r ki 011304-4 / Vol. 131, JANUARY 2009
共40兲
2
gi共r, 兲cos共m兲d
共44兲
gi共r, 兲sin共m兲d
共45兲
0 2
0
1 2
冕
2
gi共r, 兲d
共46兲
0
Similarly, Cin共兲 and Cout共兲 in Eqs. 共18兲 and 共19兲 may be expanded in a generalized Fourier series to yield boundary conditions for ODEs in Eqs. 共40兲–共42兲. Solutions for the Euler equations, Eqs. 共40兲 and 共41兲 are given by Tˆc,im共r兲 = Ac,irm + Bc,ir−m + f pc共r兲
共47兲
Transactions of the ASME
Downloaded 01 Mar 2010 to 130.126.32.13. Redistribution subject to ASME license or copyright; see http://www.asme.org/terms/Terms_Use.cfm
Fig. 3 Transient isotherms in three-layer annulus: „a… t = 5, „b… t = 10, „c… t = 15, and „d… steady state
Tˆs,im共r兲 = As,irm + Bs,ir−m + f ps共r兲
共48兲
where f pc共r兲 and f ps共r兲 are particular integrals, which can be evaluated by the application of method of variation of parameters or method of undetermined coefficients. The constants Ac,i, As,i, Bc,i, and Bs,i may be evaluated using boundary and interface conditions for Tˆc,im共r兲 and Tˆs,im共r兲. Once gˆi0共r兲 are evaluated, the solution for Tˆ 共r兲 is straightforward. i0
7
P
Illustrative Example and Results
A three-layer annulus 共r0 艋 r 艋 r3, 0 艋 艋 2; see Fig. 2兲 is initially at a uniform zero temperature. For time t ⬎ 0, -dependent heat flux given by q⬙共r = r3, 兲 =
冋
q 0 2共 − 兲 2 ,
0艋艋
0,
艋 艋 2
册
¯T 共r, ,t兲 = i
兺D
0pe
2 −␣ii0p t
Ri0p共i0pr兲
p=1
M
+
P
兺兺D
mpe
2 −␣iimp t
Rimp共impr兲cos共m兲
m=1 p=1
共49兲
is applied at the outer surface 共r = r3兲 while the inner surface 共r = r0兲 is maintained isothermal at zero temperature. This leads to the coefficients Ain = 0, Bin = 1, Aout = k3, Bout = 0, Cin共兲 = 0, and Cout共兲 = q⬙共r3 , 兲 in the respective boundary condition equations. There is no volumetric heat generation in any of the layers, i.e., gi共r , 兲 = 0. Journal of Heat Transfer
Parameter values used in this problem are k2 / k1 = 2, k3 / k1 = 4, ␣2 / ␣1 = 4, ␣3 / ␣1 = 9, r1 / r0 = 2, r2 / r0 = 4, and r3 / r0 = 6. These have been arbitrarily chosen and do not, in any way, simplify the solution. It should be noted that, in the results that follow, r, t, and Ti共r , , t兲 are in the units of r0, r20 / ␣1, and, q0r0 / k1, respectively. For this particular problem, the infinite series solution for the complementary transient temperature ¯Ti共r , , t兲 is truncated at p = P and m = M leading to
M
+
P
兺兺E
mpe
2 −␣iimp t
Rimp共impr兲sin共m兲
m=1 p=1
− i共r, ,t;M, P兲
共50兲
where i共r , , t ; M , P兲 is the truncation error. Since 1mp increases with increasing m and p, it is obvious that for a given M and P, JANUARY 2009, Vol. 131 / 011304-5
Downloaded 01 Mar 2010 to 130.126.32.13. Redistribution subject to ASME license or copyright; see http://www.asme.org/terms/Terms_Use.cfm
maximum truncation error occurs at t = 0. Moreover, since ¯T 共r , , t = 0兲 = −T 共r , 兲, therefore i ss,i P
兺D
i共r, ,t = 0;M, P兲 = Tss,i共r, 兲 +
0pRi0p共i0pr兲
p=1
P
M
+
兺兺D
mpRimp共impr兲cos共m兲
m=1 p=1 P
M
+
兺兺E
mpRimp共impr兲sin共m兲
共51兲
m=1 p=1
However, Tss,i共r , 兲 is also evaluated as a series solution 共see Eq. 共39兲兲 and hence, the above equation can be written as
冉
M ss
i共r, ,t = 0;M, P兲 = Tˆi0共r兲 +
兺 Tˆ
c,im共r兲cos共m兲
m=1 M ss
+
兺 Tˆ
s,im共r兲sin共m兲
+ ss,i共r, ,M ss兲
m=1
冊
P
+
兺D
0pRi0p共i0pr兲
p=1 M
+
P
兺兺D
mpRimp共impr兲cos共m兲
m=1 p=1 M
+
P
兺兺E
mpRimp共impr兲sin共m兲
共52兲
m=1 p=1
A good estimate of i共r , , t = 0 ; M , P兲 may be obtained only if ss,i共r , ; M ss兲 Ⰶ i共r , , t = 0 ; M , P兲. The above requirement may be fulfilled by including, not surprisingly, a large number of terms in the steady-state series solution so as to minimize the steadystate truncation error. The maximum difference between the steady-state temperatures obtained with M ss = 45 and M ss = 50 is found to be of the order of 10−5, therefore, this series is truncated at M ss = 50. The maximum percent error defined as max共共r , , t = 0兲兲 / 共Tss,max − Tss,min兲 is evaluated for various values of M and P. For M = P = 6, 8, and 10, the error is 2.25%, 1.69%, and 1.36%, respectively. Isotherms in the three-layer annulus are shown for different t values in Fig. 3. At any time t, maximum and minimum temperatures are observed on the outer edge of the annulus 共r = r3兲 at angular locations = / 2 and 3 / 2, respectively. Temperature kinks 共or discontinuity in temperature slopes兲 are clearly visible in the isotherms indicating different thermal properties of the three layers. Additionally, radial temperature variations at distinct angular positions are shown in Fig. 4. At any given 共r , t兲, maximum and minimum temperatures are observed at = / 2 and 3 / 2, respectively, as expected, since the incoming heat-flux is symmetric in 0 艋 艋 and has a maxima at = / 2. The results shown in Figs. 3 and 4 are obtained with M = P = 10.
8
Conclusions
In this paper, an analytical solution to the asymmetric transient heat conduction in a layered annulus is presented. Each layer can have spatially varying but time-independent volumetric heat source. Inhomogeneous boundary condition of the first, second, or the third kind can be applied in the radial direction. The proposed solution is also applicable to the layered-structures with inner radius r0 = 0. It is noted that the solution of the multilayer two-dimensional heat conduction problem in polar coordinates is not analogous to 011304-6 / Vol. 131, JANUARY 2009
Fig. 4 Transient temperature variation in the radial direction at „a… = 0, „b… = / 2, and „c… = 3 / 2
the corresponding problem in multidimensional Cartesian coordinates 共or 2D cylindrical r – z coordinates兲. In polar coordinates, dependence of the eigenvalues in the transverse direction on those in the other direction is not explicit. The absence of explicit dependence leads to a complete solution, which does not have imaginary transverse eigenvalues. Numerical evaluation of the Transactions of the ASME
Downloaded 01 Mar 2010 to 130.126.32.13. Redistribution subject to ASME license or copyright; see http://www.asme.org/terms/Terms_Use.cfm
double-series solution shows that a reasonable number of terms are sufficient to obtain results with acceptable errors for engineering applications.
References 关1兴 Mulholland, G. P., and Cobble, M. H., 1972, “Diffusion Through Composite Media,” Int. J. Heat Mass Transfer, 15, pp. 147–160. 关2兴 Mikhailov, M. D., Ozisik, M. N., and Vulchanov, N. L., 1983, “Diffusion in Composite Layers With Automatic Solution of the Eigenvalue Problem,” Int. J. Heat Mass Transfer, 26共8兲, pp. 1131–1141. 关3兴 Desai, A., Geer, J., and Sammakia, B., 2006, “Models of Steady Heat Conduction in Multiple Cylindrical Domains,” ASME J. Electron. Packag., 128, pp. 10–17. 关4兴 Bulavin, P. E., and Kascheev, V. M., 1965, “Solution of the NonHomogeneous Heat Conduction Equation for Multilayered Bodies,” Int. Chem. Eng., 1, pp. 112–115. 关5兴 Lu, X., Tervola, P., and Viljanen, M., 2006, “Transient Analytical Solution to Heat Conduction in Composite Circular Cylinder,” Int. J. Heat Mass Transfer, 49, pp. 341–348. 关6兴 Shupikov, A. N., Smetankina, N. V., and Svet, Y. V., 2007, “Nonstationary Heat Conduction in Complex-Shape Laminated Plates,” ASME J. Heat Transfer, 129, pp. 335–341. 关7兴 Geer, J., Desai, A., and Sammakia, B., 2007, “Heat Conduction in Multilayered Rectangular Domains,” ASME J. Electron. Packag., 129, pp. 440–451. 关8兴 Huang, S. C., and Chang, Y. P., 1980, “Heat Conduction in Unsteady, Periodic and Steady States in Laminated Composites,” ASME J. Heat Transfer, 102, pp. 742–748. 关9兴 Siegel, R., 1999, “Transient Thermal Analysis of Parallel Translucent Layers by Using Green’s Functions,” J. Thermophys. Heat Transfer, 13共1兲, pp. 10–17. 关10兴 Haji-Sheikh, A., and Beck, J. V., 1990, “Green’s Function Partitioning in Galerkin-Base Integral Solution of the Diffusion Equation,” ASME J. Heat Transfer, 112, pp. 28–34.
Journal of Heat Transfer
关11兴 Yener, Y., and Ozisik, M. N., 1974, “On the Solution of Unsteady Heat Conduction in Multi-Region Finite Media With Time-Dependent Heat Transfer Coefficient,” Proceedings of the Fifth International Heat Transfer Conference, JSME, Tokyo, Japan, Sep. 3–7, pp. 188–192. 关12兴 Salt, H., 1983, “Transient Heat Conduction in a Two-Dimensional Composite Slab. I. Theoretical Development of Temperatures Modes,” Int. J. Heat Mass Transfer, 26共11兲, pp. 1611–1616. 关13兴 Salt, H., 1983, “Transient Heat Conduction in a Two-Dimensional Composite Slab. II. Physical Interpretation of Temperatures Modes,” Int. J. Heat Mass Transfer, 26共11兲, pp. 1617–1623. 关14兴 Mikhailov, M. D., and Ozisik, M. N., 1986, “Transient Conduction in a ThreeDimensional Composite Slab,” Int. J. Heat Mass Transfer, 29共2兲, pp. 340–342. 关15兴 Haji-Sheikh, A., and Beck, J. V., 2002, “Temperature Solution in MultiDimensional Multi-Layer Bodies,” Int. J. Heat Mass Transfer, 45共9兲, pp. 1865–1877. 关16兴 de Monte, F., 2003, “Unsteady heat Conduction in Two-Dimensional Two Slab-Shaped Regions. Exact Closed-Form Solution and Results,” Int. J. Heat Mass Transfer, 46共8兲, pp. 1455–1469. 关17兴 Yeh, F.-B., 2007, “Prediction of the Transient and Steady Temperature Distributions in a Two-Layer Composite Slab in Contact With a Plasma: Exact Closed-Form Solutions,” J. Phys. D: Appl. Phys., 40, pp. 3633–3643. 关18兴 Singh, S., Jain, P. K., and Rizwan-uddin, 2008, “Analytical Solution to Transient Heat Conduction in Polar Coordinates With Multiple Layers in Radial Direction,” Int. J. Therm. Sci., 47, pp. 261–273. 关19兴 Yuen, W. Y. D., 1994, “Transient Temperature Distribution in a Multilayer Medium Subject to Radiative Surface Cooling,” Appl. Math. Model., 18, pp. 93–100. 关20兴 Singh, S., Jain, P. K., and Rizwan-uddin, 2007, “Finite Integral Transform Technique to Solve Asymmetric Heat Conduction in a Multilayer Annulus With Time-Dependent Boundary Conditions,” submitted. 关21兴 Ozisik, M. N., 1993, Heat Conduction, 2nd ed., Wiley and Sons, New York. 关22兴 de Monte, F., 2002, “An Analytic Approach to the Unsteady Heat Conduction Processes in One-Dimensional Composite Media,” Int. J. Heat Mass Transfer, 45共6兲, pp. 1333–1343.
JANUARY 2009, Vol. 131 / 011304-7
Downloaded 01 Mar 2010 to 130.126.32.13. Redistribution subject to ASME license or copyright; see http://www.asme.org/terms/Terms_Use.cfm
Analytical Solution to Transient Asymmetric Heat Conduction in a Multilayer Annulus In this paper, we present an analytical double-series solution for the time-dependent asymmetric heat conduction in a multilayer annulus. In general, analytical solutions in multidimensional Cartesian or cylindrical 共r , z兲 coordinates suffer from existence of imaginary eigenvalues and thus may lead to numerical difficulties in computing analytical solution. In contrast, the proposed analytical solution in polar coordinates (2D cylindrical) is “free” from such imaginary eigenvalues. Real eigenvalues are obtained by virtue of precluded explicit dependence of transverse (radial) eigenvalues on those in the other direction. 关DOI: 10.1115/1.2977553兴 Keywords: heat conduction, layered annulus, analytical method
1
Introduction
In modern engineering applications, multilayer components are extensively used due to the added advantage of combining physical, mechanical, and thermal properties of different materials. Many of these applications require a detailed knowledge of transient temperature and heat-flux distribution within the component layers. Both analytical and numerical techniques may be used to solve such problems. Nonetheless, numerical solutions are preferred and prevalent in practice, due to either unavailability or higher mathematical complexity of the corresponding exact solutions. Rather limited use of analytical solutions should not diminish their merit over numerical ones; since exact solutions, if available, provide an insight into the governing physics of the problem, which is typically missing in any numerical solution. Moreover, analyzing closed-form solutions to obtain optimal design options for any particular application of interest is relatively simpler. In addition, exact solutions find their applications in validating and comparing various numerical algorithms to help improve computational efficiency of computer codes that currently rely on numerical techniques. Although multilayer heat conduction problems have been studied in great detail and various solution methods—including orthogonal and quasiorthogonal expansion technique, Laplace transform method, Green’s function approach, finite integral transform technique 关1–11兴—are readily available; there is a continued need to develop and explore novel methods to solve problems for which exact solutions still do not exist. One such problem is to determine exact unsteady temperature distribution in polar coordinates 共r , 兲 with multiple layers in the radial direction. Salt 关12,13兴 addressed time-dependent heat conduction problem by orthogonal expansion technique, in a two-dimensional composite slab 共Cartesian geometry兲 with no internal heat source, subjected to homogenous boundary conditions. Later, Mikhailov and Ozisik 关14兴 solved the 3D transient conduction problem in a Cartesian nonhomogenous finite medium. More recently, Haji-Sheikh and Beck 关15兴 applied Green’s function approach to develop transient temperature solutions in a 3D Cartesian two-layer orthotropic medium including the effects of contact resistance. Lu et al. 关5兴 developed a novel method by combining Laplace transform 1 Corresponding author. Contributed by the Heat Transfer Division of ASME for publication in the JOURNAL OF HEAT TRANSFER. Manuscript received February 12, 2008; final manuscript received April 21, 2008; published online October 20, 2008. Review conducted by Peter Vadasz.
Journal of Heat Transfer
method and separation of variables method to solve multidimensional transient heat conduction problem in a rectangular and cylindrical multilayer slab with time-dependent periodic boundary condition. The treatment in the cylindrical coordinates is, however, restricted to the r – z coordinates. Eigenfunction expansion method is applied by de Monte 关16兴 to solve the unsteady heat conduction problem in a two-dimensional two-layer isotropic slab subjected to homogenous boundary conditions. Feng-Bin Yeh 关17兴 applied the method of separation of variables to solve plasma heating of a one-dimensional two-layer composite slab with layers in imperfect thermal contact. The brief review of relevant literature is by no means exhaustive. However, a literature survey showed that the analytical solution for unsteady temperature distribution in a multilayer annular geometry has not yet been developed. Recently, an exact solution based on the separation of variables method is developed by Singh et al. 关18兴 for multilayer heat conduction in polar coordinates. However, that exact solution is applicable only to domains with pie slice geometry 共 ⬍ 2, where is the angle subtended by the layers兲. Numerous applications involving multilayer cylindrical geometry require evaluation of temperature distribution in complete disk-type 共i.e., = 2兲 layers. One typical example is a nuclear fuel rod, which consists of concentric layers of different materials and often subjected to asymmetric boundary conditions. Moreover, several other applications including multilayer insulation materials, double heat-flux conductimeter, typical laser absorption calorimetry experiments, cryogenic systems, and other cylindrical building structures would benefit from such analytical solutions. This paper extends the solution approach developed by Singh et al. 关18兴 for such applications and presents an analytical double-series solution for the time-dependent asymmetric heat conduction in a multilayer annulus. Solution is valid for any combinations of time-independent, inhomogeneous boundary conditions at the inner and outer radii of the domain. The results for an illustrative problem involving a three-layer annulus subjected to asymmetric heat-flux are also presented.
2
Mathematical Formulation
Consider an n-layer annulus 共r0 艋 r 艋 rn兲, as shown schematically in Fig. 1. All the layers are assumed to be isotropic in thermal properties and are in perfect thermal contact. Let ki and ␣i be the temperature independent thermal conductivity and thermal diffusivity of the ith layer. At t = 0, each ith layer is at a specified temperature f i共r , 兲 and time-independent heat sources gi共r , 兲 are switched on for t ⬎ 0. Both the inner 共i = 1, r = r0兲 as well as the
Copyright © 2009 by ASME
JANUARY 2009, Vol. 131 / 011304-1
Downloaded 01 Mar 2010 to 130.126.32.13. Redistribution subject to ASME license or copyright; see http://www.asme.org/terms/Terms_Use.cfm
ki
Ti Ti−1 共ri−1, ,t兲 = ki−1 共ri−1, ,t兲 r r
共7兲
Initial condition 共for ri−1 艋 r 艋 ri and 0 艋 艋 2, where i = 1 , 2 , . . . , n兲: Ti共r, ,t = 0兲 = f i共r, 兲
共8兲
Boundary conditions either of the first, second, or third kind may be imposed at r = r0 and r = rn by choosing the appropriate coefficients in Eqs. 共2兲 and 共3兲. However, the case in which Bin and Bout are simultaneously zero is not considered. In addition, asymmetric boundary conditions can be applied by choosing -dependent Cin and Cout. Furthermore, multiple layers with zero inner radius 共r0 = 0兲 can be simulated by assigning zero values to constants Bin and Cin in Eq. 共2兲. It should be noted that the formulation presented in this paper only applies toward timeindependent boundary conditions and/or source terms due to the limitation of the separation of variables method. This solution methodology cannot be extended to include the effects of timedependence in boundary conditions and/or sources. Such problems can be solved analytically using the finite integral transform technique 关20,21兴. Fig. 1 Schematic of an n-layer annulus. ith layer has an inner and outer radii equal to ri−1 and ri, respectively.
outer 共i = n, r = rn兲 surfaces of the annulus may be subjected to any combination of temperature and heat-flux boundary conditions. 共Since perfect thermal contact between the adjacent layers is seldom observed in real materials, dealing with imperfect contact would require explicit modeling of the thermal resistance at the layer interfaces 关11,17,19兴. For such cases, the temperature at the contact interfaces will not be continuous.兲 Under these assumptions, the governing heat conduction equation along with the boundary and initial conditions are as follows. Governing equation 共for ri−1 艋 r 艋 ri, 0 艋 艋 2, and t ⬎ 0, where i = 1 , 2 , . . . , n兲:
冉
冊
3
Solution Methodology
In order to apply the separation of variables method, which is only applicable to homogenous problems, the nonhomogenous problem has to be split 关21兴 into: 共1兲 homogenous transient problem, and 共2兲 nonhomogenous steady-state problem. This is accomplished by splitting Ti共r , , t兲 in the governing Eqs. 共1兲–共8兲 as ¯T 共r , , t兲 + T 共r , 兲, where ¯T 共r , , t兲 is the “complementary” i
ss,i
3.1 Homogenous Transient Problem. Homogenized complementary transient equations corresponding to Eqs. 共1兲–共8兲 are as follows. Governing equation 共for ri−1 艋 r 艋 ri, 0 艋 艋 2 and t ⬎ 0, where i = 1 , 2 , . . . , n兲:
冉
1 Ti 1 2T i gi共r, 兲 1 Ti 共r, ,t兲 = r 共r, ,t兲 + 2 2 共r, ,t兲 + ␣i t r r r r ki 共1兲 Boundary conditions: •
•
Tn 共rn, ,t兲 + BoutTn共rn, ,t兲 = Cout共兲 Aout r •
•
Inner surface of the first layer 共for 0 艋 艋 2 and t ⬎ 0兲 Ain
共2兲 •
•
¯T1 共r0, ,t兲 + Bin¯T1共r0, ,t兲 = 0 r
共10兲
Outer surface of the nth layer 共for 0 艋 艋 2 and t ⬎ 0兲 Aout
共3兲
Periodic boundary conditions 共for ri−1 艋 r 艋 ri and t ⬎ 0, where i = 1 , 2 , . . . , n兲
共9兲
Boundary conditions: •
Outer surface of the nth layer 共for 0 艋 艋 2 and t ⬎ 0兲
冊
1 ¯Ti 1 2¯Ti 1 ¯Ti 共r, ,t兲 = r 共r, ,t兲 + 2 2 共r, ,t兲 ␣i t r r r r
Inner surface of the first layer 共for 0 艋 艋 2 and t ⬎ 0兲
T1 共r0, ,t兲 + BinT1共r0, ,t兲 = Cin共兲 Ain r
i
transient part and Tss,i共r , 兲 is the steady-state part of the solution.
¯Tn 共rn, ,t兲 + Bout¯Tn共rn, ,t兲 = 0 r
共11兲
Periodic boundary conditions 共for ri−1 艋 r 艋 ri and t ⬎ 0, where i = 1 , 2 , . . . , n兲
Ti共r, = 0,t兲 = Ti共r, = 2,t兲
共4兲
¯T 共r, = 0,t兲 = ¯T 共r, = 2,t兲 i i
共12兲
Ti Ti 共r, = 0,t兲 = 共r, = 2,t兲
共5兲
¯Ti ¯Ti 共r, = 0,t兲 = 共r, = 2,t兲
共13兲
Interface of the 共i − 1兲st and the ith layer 共for 0 艋 艋 2 and t ⬎ 0, where i = 2 , . . . , n兲 Ti共ri−1, ,t兲 = Ti−1共ri−1, ,t兲
011304-2 / Vol. 131, JANUARY 2009
共6兲
•
Interface of the 共i − 1兲st and the ith layer 共for 0 艋 艋 2 and t ⬎ 0 where i = 2 , . . . , n兲 ¯T 共r , ,t兲 = ¯T 共r , ,t兲 i i−1 i−1 i−1
共14兲
Transactions of the ASME
Downloaded 01 Mar 2010 to 130.126.32.13. Redistribution subject to ASME license or copyright; see http://www.asme.org/terms/Terms_Use.cfm
ki
¯Ti ¯Ti−1 共ri−1, ,t兲 = ki−1 共ri−1, ,t兲 r r
⬁
¯T 共r, ,t兲 = i
共15兲
冊
1 Tss,i 1 2Tss,i gi共r, 兲 共r, 兲 + 2 r 共r, 兲 + = 0, r r r r 2 ki ri−1 艋 r 艋 ri, 1 艋 i 艋 n
+
共17兲
Tss,1 共r0, 兲 + BinTss,1共r0, 兲 = Cin共兲 r
Tss,n Aout 共rn, 兲 + BoutTss,n共rn, 兲 = Cout共兲 r
共19兲
Tss,i Tss,i 共r, = 0兲 = 共r, = 2兲
共21兲
Tss,i Tss,i−1 共ri−1, 兲 = ki−1 共ri−1, 兲 r r
共22兲 共23兲
Solution to the Homogenous Transient Problem 4.1
¯T 共r, ,t兲 = R 共r兲⌰共兲⌫ 共t兲 i i i
冊冉
冊
1 d dRi共r兲 m r + i2 − 2 Ri共r兲 = 0 dr r dr r 2
in ri−1 艋 r 艋 ri, 共25兲
where i = 1,2, . . . ,n d2⌰共兲 + m2⌰共兲 = 0 d2 1 d⌫i共t兲 + i2⌫i共t兲 = 0 ␣i dt
in 0 艋 艋 2
Rimp共impr兲sin共m兲
共28兲
共29兲
i = 1,2, . . . ,n
The radial 共transverse兲 eigenfunction, Rimp共impr兲 in Eq. 共28兲 is
for t ⬎ 0,
兺
ki ␣i
共30兲
冕
ri
rRimp共impr兲Rimq共imqr兲dr =
ri−1
冋
0
if p ⫽ q
Nrmp if p = q
册
where Jm and Y m are Bessel functions of the first and second kind of order m, respectively. Nrmp is normalization integral in the r-direction. For the angular eigenfunctions ⌰m共兲—formed via combination of constant sin共m兲 and cos共m兲—the standard orthogonality condition is valid 关21兴. 4.3 Radial (Transverse) Eigencondition. Application of the boundary conditions 共Eqs. 共10兲 and 共11兲兲 and interface conditions 共Eqs. 共14兲 and 共15兲兲 to the transverse eigenfunction Rimp共impr兲 yields a 共2n ⫻ 2n兲 matrix for each integer value of m. Transverse eigencondition is obtained by setting the determinant of this matrix equal to zero. Roots of which, in turn, yield the infinite number of eigenvalues 1mp corresponding to the first layer for each integer value of m. 4.4 Determination of Coefficients aimp and bimp. Coefficients aimp and bimp in the radial eigenfunction Rimp共impr兲 共see Eq. 共30兲兲 are evaluated from the following recurrence relationship, obtained from the ith interface condition, valid for i 苸 共1 , n − 1兲:
冉 冊冉 bi+1,mp
=
Jm共i+1,mpri兲
冊 冊冉 冊
Y m共i+1,mpri兲
−1
⬘ 共i+1,mpri兲 ki+1Y m⬘ 共i+1,mpri兲 ki+1Jm
⫻
冉
Jm共impri兲
Y m共impri兲
⬘ 共impri兲 kiY m⬘ 共impri兲 k iJ m
aimp
bimp
共32兲
where b1mp = −共C1in / C2in兲a1mp, and a1mp is arbitrary. 4.5 Determination of Coefficients D0p, Dmp, and Emp. Coefficients D0p, Dmp, and Emp in Eq. 共28兲 are evaluated by applying the initial condition 共Eq. 共16兲兲 and then making use of the orthogonality conditions in the radial and angular directions, as follows n
D0p =
兺
1 ki 2Nr0p i=1 ␣i
冕冕 2
0
ri
¯ 共r, ,t = 0兲drd rRi0p共i0pr兲T i
ri−1
where i = 1,2, . . . ,n
共33兲
are constants of separation.
4.2 General Solution. In view of the ODEs listed above, a general solution for Eq. 共9兲 may be realized as Journal of Heat Transfer
n
共26兲
共27兲 where
2 −␣iimp t
where continuity of the heat-flux at the layer interfaces requires the following relationship between the ith imp and 1mp to hold 关15,16,21,22兴,
共24兲
in Eq. 共9兲 and applying separation of variables yield the following ordinary differential equations 共ODEs兲.
2i
mpe
m=1 p=1
ai+1,mp
Separation of Variables. Substituting the product form
冉
兺兺E
共31兲
共20兲
ki
Rimp共impr兲cos共m兲
⬁
⬁
+
i=1
Tss,i共r, = 0兲 = Tss,i共r, = 2兲
Tss,i共ri−1, 兲 = Tss,i−1共ri−1, 兲
4
共18兲
Interface of the 共i − 1兲st and the ith layer 共for 0 艋 艋 2, where i = 2 , . . . , n兲
•
2 −␣iimp t
and the corresponding orthogonality condition is 关18兴
Periodic boundary conditions 共for ri−1 艋 r 艋 ri, where i = 1 , 2 , . . . , n兲
•
mpe
Rimp共impr兲 = aimpJm共impr兲 + bimpY m共impr兲
Outer surface of the nth layer 共for 0 艋 艋 2兲
•
兺兺D
imp = 1mp冑␣1/␣i
Inner surface of the first layer 共for 0 艋 艋 2兲 Ain
Ri0p共i0pr兲
m=1 p=1
Boundary conditions: •
2 −␣ii0p t
⬁
⬁
共16兲
3.2 Inhomogeneous Steady-State Problem. Inhomogeneous steady-state equations corresponding to Eqs. 共1兲–共8兲 are as follows. Governing equation 共for ri−1 艋 r 艋 ri and 0 艋 艋 2, where i = 1 , 2 , . . . , n兲:
冉
0pe
p=1
Initial condition 共for ri−1 艋 r 艋 ri and 0 艋 艋 2, where i = 1 , 2 , . . . , n兲: ¯T 共r, ,t = 0兲 = f 共r, 兲 − T 共r, 兲 i i ss,i
兺D
n
兺
1 ki Dmp = Nrmp i=1 ␣i
冕冕 2
0
¯ 共r, ,t = 0兲drd ⫻T i
ri
rRimp共impr兲cos共m兲
ri−1
共34兲
JANUARY 2009, Vol. 131 / 011304-3
Downloaded 01 Mar 2010 to 130.126.32.13. Redistribution subject to ASME license or copyright; see http://www.asme.org/terms/Terms_Use.cfm
n
Emp =
兺
1 ki Nrmp i=1 ␣i
冕冕 2
0
ri
¯ 共r, ,t = 0兲drd rRimp共impr兲sin共m兲T i
ri−1
共35兲
5
Absence of imaginary Radial Eigenvalues
In general, transverse eigenvalues for multilayer timedependent heat conduction problems in Cartesian coordinates may be imaginary 关16兴. Same is true for 2D 共r , z兲 cylindrical coordinates. The eigenvalues are imaginary due to the explicit dependence of the transverse eigenvalues on those in the remaining direction共s兲. For example, in a two-layer homogenous heat conduction problem with layers in the x-direction, the general solution in the ith layer is 关15,16兴 ⬁
⬁
Ti共x,y,t兲 =
兺兺Z
mpe
2 2 −␣i共imp +m 兲t
Ximp共impx兲Y m共my兲 共36兲
m=1 p=1
where imp and m are eigenvalues in the x- and y-directions, respectively; and Ximp 共impx兲 and Y m共my兲 are the corresponding eigenfunctions. For heat-flux to be continuous at the interface ∀t 2 2 2 2 ␣1共1mp + m 兲 = ␣2共2mp + m 兲
which implies
2mp =
冑
共37兲
冉 冊
␣1 2 ␣1 2 + − 1 m ␣2 1mp ␣2
共38兲
Clearly, the relationship above may result in either real or imaginary transverse eigenvalues 关13,15,16兴. However, in the present case, similar considerations led to Eq. 共29兲, which is similar to what has been established for 1D multilayer time-dependent problems and eliminates the possibility of imaginary eigenvalues. Moreover, physical considerations dictate that eigenvalues should be real 共both in the present and the corresponding 1D case兲 because imaginary eigenvalues will result in exponentially growing temperatures. In contrast, the solution given in Eq. 共36兲 can have a physically realizable solution even in the case where eigenvalues in one of the directions are imaginary. It should be noted that, though there is no explicit dependence between radial and angular eigenvalues, the order of the Bessel functions constituting radial eigenfunctions is determined by the angular eigenvalues. Hence, the radial eigenvalues implicitly depend on the angular eigenvalues. Moreover, unlike in Cartesian coordinates, this implicit dependence does not vanish even if ␣1 = ␣i共i ⫽ 1兲. In fact, it exists even for single-layer problems.
Fig. 2 Asymmetric heat conduction in a three-layer annulus. Each layer has a different thermal conductivity „ki… and thermal diffusivity „␣i…. The lower-half of the annulus „ Ï Ï 2… is kept insulated, while the upper-half „0 Ï Ï … is subjected to a -dependent incoming heat-flux.
冉
冊
1 d dTˆs,im共r兲 m2 gˆs,im共r兲 r − 2 Tˆs,im共r兲 + =0 dr r dr r ki
冉
冊
1 d dTˆi0共r兲 gˆi0共r兲 r + =0 dr r dr ki
共41兲
共42兲
where the source term gi共r , 兲 is expanded in a generalized Fourier series as ⬁
gi共r, 兲 = gˆi0共r兲 +
兺
⬁
gˆc,im共r兲cos共m兲 +
m=1
兺 gˆ
s,im共r兲sin共m兲
m=1
共43兲 and
6 Solution to the Inhomogeneous Steady-State Problem The inhomogeneous steady-state problem is solved using the eigenfunction expansion method. The steady-state temperature distribution, governed by Eq. 共17兲, may be written as a generalized Fourier series in terms of angular eigenfunctions, ⬁
Tss,i共r, 兲 = Tˆi0共r兲 +
兺 Tˆ
s,im共r兲sin共m兲,
ri−1 艋 r 艋 ri,
1艋i艋n
m=1
共39兲 Substituting Eq. 共39兲 in Eq. 共17兲 leads to the following ODEs:
冉
冕
1
冕
gˆi0共r兲 ⬅
⬁
兺 Tˆ
gˆs,im共r兲 ⬅
1
c,im共r兲cos共m兲
m=1
+
gˆc,im共r兲 ⬅
冊
1 d dTˆc,im共r兲 m2 gˆc,im共r兲 r − 2 Tˆc,im共r兲 + =0 dr r dr r ki 011304-4 / Vol. 131, JANUARY 2009
共40兲
2
gi共r, 兲cos共m兲d
共44兲
gi共r, 兲sin共m兲d
共45兲
0 2
0
1 2
冕
2
gi共r, 兲d
共46兲
0
Similarly, Cin共兲 and Cout共兲 in Eqs. 共18兲 and 共19兲 may be expanded in a generalized Fourier series to yield boundary conditions for ODEs in Eqs. 共40兲–共42兲. Solutions for the Euler equations, Eqs. 共40兲 and 共41兲 are given by Tˆc,im共r兲 = Ac,irm + Bc,ir−m + f pc共r兲
共47兲
Transactions of the ASME
Downloaded 01 Mar 2010 to 130.126.32.13. Redistribution subject to ASME license or copyright; see http://www.asme.org/terms/Terms_Use.cfm
Fig. 3 Transient isotherms in three-layer annulus: „a… t = 5, „b… t = 10, „c… t = 15, and „d… steady state
Tˆs,im共r兲 = As,irm + Bs,ir−m + f ps共r兲
共48兲
where f pc共r兲 and f ps共r兲 are particular integrals, which can be evaluated by the application of method of variation of parameters or method of undetermined coefficients. The constants Ac,i, As,i, Bc,i, and Bs,i may be evaluated using boundary and interface conditions for Tˆc,im共r兲 and Tˆs,im共r兲. Once gˆi0共r兲 are evaluated, the solution for Tˆ 共r兲 is straightforward. i0
7
P
Illustrative Example and Results
A three-layer annulus 共r0 艋 r 艋 r3, 0 艋 艋 2; see Fig. 2兲 is initially at a uniform zero temperature. For time t ⬎ 0, -dependent heat flux given by q⬙共r = r3, 兲 =
冋
q 0 2共 − 兲 2 ,
0艋艋
0,
艋 艋 2
册
¯T 共r, ,t兲 = i
兺D
0pe
2 −␣ii0p t
Ri0p共i0pr兲
p=1
M
+
P
兺兺D
mpe
2 −␣iimp t
Rimp共impr兲cos共m兲
m=1 p=1
共49兲
is applied at the outer surface 共r = r3兲 while the inner surface 共r = r0兲 is maintained isothermal at zero temperature. This leads to the coefficients Ain = 0, Bin = 1, Aout = k3, Bout = 0, Cin共兲 = 0, and Cout共兲 = q⬙共r3 , 兲 in the respective boundary condition equations. There is no volumetric heat generation in any of the layers, i.e., gi共r , 兲 = 0. Journal of Heat Transfer
Parameter values used in this problem are k2 / k1 = 2, k3 / k1 = 4, ␣2 / ␣1 = 4, ␣3 / ␣1 = 9, r1 / r0 = 2, r2 / r0 = 4, and r3 / r0 = 6. These have been arbitrarily chosen and do not, in any way, simplify the solution. It should be noted that, in the results that follow, r, t, and Ti共r , , t兲 are in the units of r0, r20 / ␣1, and, q0r0 / k1, respectively. For this particular problem, the infinite series solution for the complementary transient temperature ¯Ti共r , , t兲 is truncated at p = P and m = M leading to
M
+
P
兺兺E
mpe
2 −␣iimp t
Rimp共impr兲sin共m兲
m=1 p=1
− i共r, ,t;M, P兲
共50兲
where i共r , , t ; M , P兲 is the truncation error. Since 1mp increases with increasing m and p, it is obvious that for a given M and P, JANUARY 2009, Vol. 131 / 011304-5
Downloaded 01 Mar 2010 to 130.126.32.13. Redistribution subject to ASME license or copyright; see http://www.asme.org/terms/Terms_Use.cfm
maximum truncation error occurs at t = 0. Moreover, since ¯T 共r , , t = 0兲 = −T 共r , 兲, therefore i ss,i P
兺D
i共r, ,t = 0;M, P兲 = Tss,i共r, 兲 +
0pRi0p共i0pr兲
p=1
P
M
+
兺兺D
mpRimp共impr兲cos共m兲
m=1 p=1 P
M
+
兺兺E
mpRimp共impr兲sin共m兲
共51兲
m=1 p=1
However, Tss,i共r , 兲 is also evaluated as a series solution 共see Eq. 共39兲兲 and hence, the above equation can be written as
冉
M ss
i共r, ,t = 0;M, P兲 = Tˆi0共r兲 +
兺 Tˆ
c,im共r兲cos共m兲
m=1 M ss
+
兺 Tˆ
s,im共r兲sin共m兲
+ ss,i共r, ,M ss兲
m=1
冊
P
+
兺D
0pRi0p共i0pr兲
p=1 M
+
P
兺兺D
mpRimp共impr兲cos共m兲
m=1 p=1 M
+
P
兺兺E
mpRimp共impr兲sin共m兲
共52兲
m=1 p=1
A good estimate of i共r , , t = 0 ; M , P兲 may be obtained only if ss,i共r , ; M ss兲 Ⰶ i共r , , t = 0 ; M , P兲. The above requirement may be fulfilled by including, not surprisingly, a large number of terms in the steady-state series solution so as to minimize the steadystate truncation error. The maximum difference between the steady-state temperatures obtained with M ss = 45 and M ss = 50 is found to be of the order of 10−5, therefore, this series is truncated at M ss = 50. The maximum percent error defined as max共共r , , t = 0兲兲 / 共Tss,max − Tss,min兲 is evaluated for various values of M and P. For M = P = 6, 8, and 10, the error is 2.25%, 1.69%, and 1.36%, respectively. Isotherms in the three-layer annulus are shown for different t values in Fig. 3. At any time t, maximum and minimum temperatures are observed on the outer edge of the annulus 共r = r3兲 at angular locations = / 2 and 3 / 2, respectively. Temperature kinks 共or discontinuity in temperature slopes兲 are clearly visible in the isotherms indicating different thermal properties of the three layers. Additionally, radial temperature variations at distinct angular positions are shown in Fig. 4. At any given 共r , t兲, maximum and minimum temperatures are observed at = / 2 and 3 / 2, respectively, as expected, since the incoming heat-flux is symmetric in 0 艋 艋 and has a maxima at = / 2. The results shown in Figs. 3 and 4 are obtained with M = P = 10.
8
Conclusions
In this paper, an analytical solution to the asymmetric transient heat conduction in a layered annulus is presented. Each layer can have spatially varying but time-independent volumetric heat source. Inhomogeneous boundary condition of the first, second, or the third kind can be applied in the radial direction. The proposed solution is also applicable to the layered-structures with inner radius r0 = 0. It is noted that the solution of the multilayer two-dimensional heat conduction problem in polar coordinates is not analogous to 011304-6 / Vol. 131, JANUARY 2009
Fig. 4 Transient temperature variation in the radial direction at „a… = 0, „b… = / 2, and „c… = 3 / 2
the corresponding problem in multidimensional Cartesian coordinates 共or 2D cylindrical r – z coordinates兲. In polar coordinates, dependence of the eigenvalues in the transverse direction on those in the other direction is not explicit. The absence of explicit dependence leads to a complete solution, which does not have imaginary transverse eigenvalues. Numerical evaluation of the Transactions of the ASME
Downloaded 01 Mar 2010 to 130.126.32.13. Redistribution subject to ASME license or copyright; see http://www.asme.org/terms/Terms_Use.cfm
double-series solution shows that a reasonable number of terms are sufficient to obtain results with acceptable errors for engineering applications.
References 关1兴 Mulholland, G. P., and Cobble, M. H., 1972, “Diffusion Through Composite Media,” Int. J. Heat Mass Transfer, 15, pp. 147–160. 关2兴 Mikhailov, M. D., Ozisik, M. N., and Vulchanov, N. L., 1983, “Diffusion in Composite Layers With Automatic Solution of the Eigenvalue Problem,” Int. J. Heat Mass Transfer, 26共8兲, pp. 1131–1141. 关3兴 Desai, A., Geer, J., and Sammakia, B., 2006, “Models of Steady Heat Conduction in Multiple Cylindrical Domains,” ASME J. Electron. Packag., 128, pp. 10–17. 关4兴 Bulavin, P. E., and Kascheev, V. M., 1965, “Solution of the NonHomogeneous Heat Conduction Equation for Multilayered Bodies,” Int. Chem. Eng., 1, pp. 112–115. 关5兴 Lu, X., Tervola, P., and Viljanen, M., 2006, “Transient Analytical Solution to Heat Conduction in Composite Circular Cylinder,” Int. J. Heat Mass Transfer, 49, pp. 341–348. 关6兴 Shupikov, A. N., Smetankina, N. V., and Svet, Y. V., 2007, “Nonstationary Heat Conduction in Complex-Shape Laminated Plates,” ASME J. Heat Transfer, 129, pp. 335–341. 关7兴 Geer, J., Desai, A., and Sammakia, B., 2007, “Heat Conduction in Multilayered Rectangular Domains,” ASME J. Electron. Packag., 129, pp. 440–451. 关8兴 Huang, S. C., and Chang, Y. P., 1980, “Heat Conduction in Unsteady, Periodic and Steady States in Laminated Composites,” ASME J. Heat Transfer, 102, pp. 742–748. 关9兴 Siegel, R., 1999, “Transient Thermal Analysis of Parallel Translucent Layers by Using Green’s Functions,” J. Thermophys. Heat Transfer, 13共1兲, pp. 10–17. 关10兴 Haji-Sheikh, A., and Beck, J. V., 1990, “Green’s Function Partitioning in Galerkin-Base Integral Solution of the Diffusion Equation,” ASME J. Heat Transfer, 112, pp. 28–34.
Journal of Heat Transfer
关11兴 Yener, Y., and Ozisik, M. N., 1974, “On the Solution of Unsteady Heat Conduction in Multi-Region Finite Media With Time-Dependent Heat Transfer Coefficient,” Proceedings of the Fifth International Heat Transfer Conference, JSME, Tokyo, Japan, Sep. 3–7, pp. 188–192. 关12兴 Salt, H., 1983, “Transient Heat Conduction in a Two-Dimensional Composite Slab. I. Theoretical Development of Temperatures Modes,” Int. J. Heat Mass Transfer, 26共11兲, pp. 1611–1616. 关13兴 Salt, H., 1983, “Transient Heat Conduction in a Two-Dimensional Composite Slab. II. Physical Interpretation of Temperatures Modes,” Int. J. Heat Mass Transfer, 26共11兲, pp. 1617–1623. 关14兴 Mikhailov, M. D., and Ozisik, M. N., 1986, “Transient Conduction in a ThreeDimensional Composite Slab,” Int. J. Heat Mass Transfer, 29共2兲, pp. 340–342. 关15兴 Haji-Sheikh, A., and Beck, J. V., 2002, “Temperature Solution in MultiDimensional Multi-Layer Bodies,” Int. J. Heat Mass Transfer, 45共9兲, pp. 1865–1877. 关16兴 de Monte, F., 2003, “Unsteady heat Conduction in Two-Dimensional Two Slab-Shaped Regions. Exact Closed-Form Solution and Results,” Int. J. Heat Mass Transfer, 46共8兲, pp. 1455–1469. 关17兴 Yeh, F.-B., 2007, “Prediction of the Transient and Steady Temperature Distributions in a Two-Layer Composite Slab in Contact With a Plasma: Exact Closed-Form Solutions,” J. Phys. D: Appl. Phys., 40, pp. 3633–3643. 关18兴 Singh, S., Jain, P. K., and Rizwan-uddin, 2008, “Analytical Solution to Transient Heat Conduction in Polar Coordinates With Multiple Layers in Radial Direction,” Int. J. Therm. Sci., 47, pp. 261–273. 关19兴 Yuen, W. Y. D., 1994, “Transient Temperature Distribution in a Multilayer Medium Subject to Radiative Surface Cooling,” Appl. Math. Model., 18, pp. 93–100. 关20兴 Singh, S., Jain, P. K., and Rizwan-uddin, 2007, “Finite Integral Transform Technique to Solve Asymmetric Heat Conduction in a Multilayer Annulus With Time-Dependent Boundary Conditions,” submitted. 关21兴 Ozisik, M. N., 1993, Heat Conduction, 2nd ed., Wiley and Sons, New York. 关22兴 de Monte, F., 2002, “An Analytic Approach to the Unsteady Heat Conduction Processes in One-Dimensional Composite Media,” Int. J. Heat Mass Transfer, 45共6兲, pp. 1333–1343.
JANUARY 2009, Vol. 131 / 011304-7
Downloaded 01 Mar 2010 to 130.126.32.13. Redistribution subject to ASME license or copyright; see http://www.asme.org/terms/Terms_Use.cfm
