GREEN’S FUNCTION FOR HEAT CONDUCTION PROBLEMS IN A MULTI-LAYERED HOLLOW CYLINDER
Stanisław Kukla, Urszula Siedlecka
Institute of Mathematics, Czestochowa University of Technology Częstochowa, Poland
stanislaw.kukla@im.pcz.pl, urszula.siedlecka@im.pcz.pl
Abstract. In this paper, derivation of the Green’s function for the heat conduction problems in a finite multi-layered hollow cylinder is presented. Formulation and solution of the prob- lem includes an arbitrary number of the cylinder layers characterized by various thermal properties. At the interfaces perfect thermal contact was assumed. The Green’s function for the three-dimensional heat conduction problems in the cylindrical coordinate has been presented in the form of a product of two other Green’s functions.
Keywords: Green’s function, heat conduction, multi-layered composite cylinder
Introduction
The Green’s function (GF) method has been widely used in the solution of heat conduction problems for basic geometries in the considered regions in the cases of classical boundary conditions. Application of the Green’s function to such prob- lems can be found in the book by Beck et al. [1], in the book by Özişik [2] and in the book by Duffy [3], as well as in several papers, for example in references [4-7].
The advantage of the Green’s function method lies in the fact that the solution of the non-homogeneous problems can be expressed in terms of the Green’s function.
In the solution to the problem the terms which influence the temperature distribu- tion can be distinguished: (1) the term which gives the contribution of the initial condition, (2) the term which represents the contribution of the energy generation and (3) the term which represents the contribution of the non-homogeneous bound- ary conditions.
In the books [1-3] and in the papers [4-7] derivations of the Green’s functions for the heat conduction equation for various boundary conditions are presented and examples of applications are given. The problems presented in the references [1-7] concern the systems characterized by constant thermomechanical quantities.
In the book by Özişik [3] and in the paper by Milligan and Kinra [4] an application of the GF method to the one-dimensional heat conduction in a layered medium is also given.
The purpose of this study is to develop an analytical solution of the three- -dimensional heat conduction problem in a hollow multilayered cylinder by using
the Green’s function method. The formulation and solution to the problem concern the finite cylinder with the boundary condition of a third kind and perfect thermal contact at the interfaces.
1. Formulation of the problem
Consider a finite, N-layered hollow cylinder as is shown in Figure 1. In order to solve the heat conduction problem in the cylinder by using the GF method, an auxiliary problem is formulated [1]. This problem lies on the derivation of the necessary Green’s functions Gij, ,i j=1, 2,...,N. The functions Gij in the cylin- drical coordinates satisfy the following differential equation
( ) ( ) ( ) ( )
[ ] [ ]
2 2 2
2 2 2 2 0
1 1
1 1 1 1
0, , , , , , 1, 2,...,
i j i j i j ij i j
i i
i i j j
G G G G G
r z t
r r r t
r r z
z L r r r r r i j N
δ ρ δ ϕ ϕ δ ζ δ τ
λ α
ϕ
− ρ −
∂ ∂ ∂ ∂ ∂
+ + + + − − − − =
∂ ∂
∂ ∂ ∂
∈ ∈ ∈ =
(1) where αi and λi are the thermal diffusivity and thermal conductivity, respectively,
0
, , , , ,
r ϕ z ρ ϕ ζ are the cylindrical coordinates, δ is the Dirac delta function.
Fig. 1. The sketch of the N-layered hollow cylinder
We assume zero initial condition and convective conditions (the third kind of boundaryconditions forthetemperature distribution)at the boundaries of the hollow cylinder. The Green’s functions Gij
(
t r, ,ϕ, ; ,zτ ρ ϕ ζ, 0,)
associated with such a heat conduction problem satisfy the following homogeneous boundary conditions* 1
1 j 1 for 0
j
G
G r r
r
λ α
+
∂
− = − =
∂ (2)
* N j for
N N j N
G
G r r
r
λ α
∞
∂
− = =
∂ (3)
L
r
0
r
N
1 i j 2 for 0
i j
G
G z
γ r γ
∂
= =
∂ (4)
1 i j 2 for
i j
G
G z L
r
ν ν
∂
= =
∂ (5)
where λ1*,λN*,γ ν are the thermal conductivities, 1, 1 α+, α∞,γ2,ν2 are the heat transfer coefficients. Moreover, the conditions of perfect thermal contact at the interfaces holds [2]
1 for , , 1, 2,..., 1
i j i j i
G =G+ r=r i j= N− (6)
1
1 for , , 1, 2,..., 1
i j i j
i i i
G G
r r i j N
r r
λ λ +
+
∂ ∂
= = = −
∂ ∂ (7)
The GF for heat conduction in a cylinder with constant thermal properties (dif- fusivity and conductivity) can be presented in the form of a product of two Green’s functions [1]. Similarly, in the case of N-layered cylinder, each function Gi j can be expressed in the form of the product
(
, , , ; , , 0,) (
, , ; , , 0) (
, ; ,)
i j i j
G t rϕ zτ ρ ϕ ζ =G% t rϕ τ ρ ϕ ⋅G t zτ ζ (8) The Green’s functions for the radial direction heat conduction G% and for one- i j dimensional heat conduction G are solutions of the differential equations
( ) ( ) ( )
2 2
2 2 2 0
1 1 1 1
, 1, 2,...,
i j i j i j i j
i i
G G G G
r t
r r r t
r r
i j N
δ ρ δ ϕ ϕ δ τ
λ α
ϕ
∂ ∂ ∂ ∂
+ + + − − − =
∂ ∂
∂ ∂
=
% % % %
(9)
( ) ( )
2 2
G G
z t
z δ ζ δ τ t
∂ ∂
+ − − =
∂ ∂ (10)
The boundary conditions for the functions and
G%i j G are obtained by substituting the product (8) in equations (2)-(7).
2. Derivation of the Green’s functions
The functions G% occurring in equation (9) we find in the form of the cosine i j series:
(
0) ( ) (
0)
0
, , ; , , 1 , ; , cos
i j i j n
n n
G t rϕ τ ρ ϕ g t rτ ρ n ϕ ϕ
κ
∞
=
=
∑
−% (11)
where κ0=2π and
κn =π for n =1, 2,.... Taking into account the series (11) in the equation (9), we obtain an equation for the radial direction Green’s functions in the form
( ) ( )
[ ]
2 2
2 2
1 1
1 1 1
,
, , , , , 1, 2,...,
i j n i j n
i i
i i j j
n g
g r t
r r r t
r r
r r r r r i j N
δ ρ δ τ
λ α
− ρ −
∂ ∂ ∂
+ − + − − =
∂ ∂
∂
∈ ∈ =
(12)
2.1. Derivation of the radial direction Green’s functions
The functions gijn
(
t r τ ρ, ; ,)
satisfy the differential equation (12), the zero initial condition and the boundary conditions which follow from conditions (2)-(3):* 1
1 j n 1 for 0
j n
g
g r r
r
λ α
+
∂
− = − =
∂ (13)
* N j n for
N N j n N
g
g r r
r
λ α
∞
∂
− = =
∂ (14)
Moreover, on the basis of equations (6)-(7), using (8) and (11), we obtain:
1 for , 1, 2,..., 1
ij n i j n i
g =g+ r=r i= N− (15)
1
1 for , 1, 2,..., 1
ij n i j n
i i i
g g
r r i N
r r
λ λ +
+
∂ ∂
= = = −
∂ ∂ (16)
The functions gijn
(
t r τ ρ, ; ,)
we present in the form of an expansion( ) ( ) ( ) ( )
1
, ; , ; imn jmn
i j n mn r
m mn
g t r t r
N
φ φ ρ
τ ρ τ
∞
=
=
∑
Θ (17)where φimn
( )
r are the eigenfunctions of the following boundary problem[ ]
2 2
2
2 2 1
1 im n n i m n 0, i , i , 1, 2,...,
r r r i N
r r
r µ r φ
−
∂ ∂
+ + − = ∈ =
∂ ∂
(18)
( ) ( )
* 1
1 0 1 0
m n
m n
d
r r
dr
λ φ α φ+
− = − (19)
( )
1( )
, 1, 2,...,im n ri i mn ri i N
φ =φ+ = (20)
( )
1 1( )
, 1, 2,...,im n i m n
i i i i
d d
r r i N
dr dr
φ φ
λ =λ+ + = (21)
( ) ( )
* Nm n
N N Nm n N
d
r r
dr
λ φ α φ∞
− = (22)
where imn2 mn2
i
γ µ
= α .
The general solution of equation (18) has the form
( )
1( )
2( )
,[
1,]
, 1, 2,...,i m n r ci Jn im nr ci Yn im nr r ri ri i N
φ µ µ
= + ∈ − = (23)
where and
n n
J Y are the Bessel functions of the first and second kind, respective- ly. Substituting the functions (23) into boundary conditions (19)-(22) we obtain a system of 2N homogeneous equations. We have written the equation system in the matrix form
AC=0 (24)
where
1 , 2 ij i j N
A ≤ ≤
=
A , 11 12 21 22 ... ,1 ,2
T
N N
c c c c c c
=
C . Thenon-zeroelements
of the matrix A are 11 12 21 22 ... ,1 ,2 T
N N
c c c c c c
=
C
( ) ( )
11 * 1 0 1 1 0
1 0 1 1
n m n n m n
mn mn
A n J r J r
r
α µ µ
µ λ µ
+
+
= − −
,
( ) ( )
12 * 1 0 1 1 0
1 0 1 1
n m n n m n
mn mn
A n Y r Y r
r
α µ µ
µ λ µ
+
+
= − −
,
( ) ( )
2 ,2i i 1 n im n i , 2 ,2i i n im n i
A J µ r A Y µ r
− = = ,
( ) ( )
2 ,2i i 1 n i 1m n i , 2 ,2i i 2 n i 1m n i
A + = −J µ+ r A + = −Y µ+ r ,
( ) ( ) ( ) ( )
2i 1,2i 1 n im n i n 1 im n i , 2i 1,2i n im n i n 1 im n i
i imn i imn
n n
A J r J r A Y r Y r
r µ µ r µ µ
µ µ
+ − = − + + = − +
( ) ( )
1 1
2 1,2 1 1 1 1
1
i i mn
i i n i m n i n i m n i
i imn i i mn
n
A J r J r
r
λ µ µ µ
λ µ µ
+ +
+ + + + +
+
= − −
,
( ) ( )
1 1
2 1,2 2 1 1 1
1 1
i i mn
i i n i m n i n i m n i
i imn i i mn
A n Y r Y r
r
λ µ
µ µ
λ µ µ
+ +
+ + + + +
+ +
= − −
for i=1, 2,...,N− , 1
( ) ( )
2N,2N 1 * n Nm n N n 1 Nm n N
N Nmn N Nmn
A n J r J r
r
α µ µ
µ λ µ
∞
− +
= + −
and
( ) ( )
2N,2N * n Nm n N n 1 Nm n N
N Nmn N Nmn
A n Y r Y r
r
α µ µ
µ λ µ
∞
+
= + −
.
The non-zero solution of the system (24) exists for values of γmn for which the determinant of the matrix A is equal to zero:
detA=0 (25)
The equation (25) with respect to the eigenvalues γmn is then solved numerically.
The eigenfunctions φi m n corresponding to the calculated eigenvalues will be fully determined if the constants c11, c12, c21, c22, ...,cN,1, cN,2, occurring in equa- tions (23) will be appointed. In determining these constants we assume cN,2=1 in equation (24) and delete the last row in the matrix A in this equation. As a result we obtain the matrix equation in the form
AC=B (26)
where the matrix A arises from the matrix A by deletion of the last row and the last column, 11 12 21 22 ... ,1
T
c c c c cN
=
C and =B
(
1) (
1)
1(
1)
1 1
0 0 0 0 ...
T
N Nmn
n Nm n N n Nm n N n Nm n N
N N mn N Nmn
Y r n Y r Y r
r
µ λ µ µ µ
λ µ µ
− − + −
− −
−
The functions φi m n
( )
r , m = 1, 2,…, given by equation (23) with 1, 2i i
c c deter- mined from equation (26) create an orthogonal set of functions for any n = 1, 2,…, i.e. the following condition holds
( ) ( )
1
' 1
0 for '
for '
i
i N r
i
i m n i m n r
i i r mn
m m
r r r dr
N m m
λ α
= −
≠
Φ Φ =
=
∑ ∫
(27)The differential equation for the function Θmn
(
t τ;)
occuring in equation (17) is obtained by substituting the series (17) into equation (12) and using the orthogonal- ity condition (27). As a result we obtain the equation in the form( )
( ) ( )
; 2 mn ;
mn mn
t
t t
t
τ γ τ δ τ
∂Θ
+ Θ = −
∂ (28)
Solving equation (28) taking into account the zero initial condition one obtains
(
;)
exp 2( ) ( )
mn tτ γmn t τ H t τ
Θ = − − − (29)
Hence, the radial direction Green’s function sounds like
( ) ( )
2( ) ( ) ( )
1
, ; , exp imn jmn
i j n mn r
m mn
g t r H t t r
N
φ φ ρ
τ ρ τ γ τ
∞
=
= −
∑
− − (30)where functions φi m n
( )
r are given by equation (23).2.2. Green’s function for one-dimensional heat conduction
The Green’s function G t z
(
, ; ,τ ζ)
occurring in equation (8) is a solution of the differential equation (10) and satisfies the zero initial condition and the boundary conditions which follow from conditions (4)-(5):1 G 2 for 0
G z
γ z γ
∂ = =
∂ (31)
1 G 2 for
G z L
z
ν ν
∂ = =
∂ (32)
This function is given in the book [3] in the form of infinite series and it is presented below in completion of the derived Green’s function for heat conduction in the finite hollow cylinder. The function is given by
( ) (
2( ) ) ( ) ( )
1
, ; , exp k k k
k kz
z
G t z t
N
τ ζ β τ ζ
∞
=
Ψ Ψ
=
∑
− − (33)where
( )
1 cos 2sink z k kz kz
ψ =γ β β + γ β (34)
( )
(
1 2 12 2 22)
1 2( )
2 12 2 22( )
1 1
cos 2 sin sin 2
2 2 4
z k
k k k k k
k
N L L z γ β γ L
γ γ γ β γ γ γ β γ β β
β
= + + − + + −
and βk are roots of equation
(
1 2 2 1)
cos(
1 1 2 2 2)
sin 0k kL k kL
γ ν −γ ν β β + γ ν β +γ ν β = (35)
Finally, the Green’s function for the heat conduction problems in the finite multi-layer hollow cylinder with boundary conditions of the third kind and perfect contact at the interfaces has the form
( )
( ) ( ) ( ) ( ) ( ) ( ) ( )
0
2 2
0
1 1 1
, , , ; , , ,
exp 1 cos
i j
imn jmn k k
mn k r z
n k m mn n n
G t r z
r z
t n
N N
ϕ τ ρ ϕ ζ
ϕ ϕ ρ ζ
γ ϑ τ ϕ ϕ
κ
∞ ∞ ∞
= = =
=
Ψ Ψ
=
∑ ∑∑
− + − −(36)
Conclusions
The Green’s function for three-dimensional heat conduction problems in a finite hollow N-layered cylinder with perfect thermal contact at interfaces has been derived. The function is necessary in solving a nonhomogenous heat conduction problem in the cylinder by using the Green’s function method. The function is pre- sented in the form of a product of two Green’s functions. To determine these func- tions two eigenproblems must be solved. Use of the package Mathematica was a great assistance in symbolic computations to derive the Green’s function. The use of the Green's function defined here will assist the study of heat flow in composite cylinders.
References
[1] Beck J.V., Cole K.D., Haji-Sheikh A., Litkouhi B., Heat Conduction Using Green’s Functions, Hemisphere, Washington DC 1992.
[2] Özişik M.N., Heat Condition, second edition, John Wiley & Sons, Inc., New York 1993.
[3] Duffy D.G., Green’s Functions with Applications, Chapman&Hall/CRC, Washington DC 2001.
[4] Milligan K.B., Kinra V.K., Elastothermodynamic damping of fiber-reinforced metal-matrix composites, Journal of Applied Mechanics 1995, 62, 441-449.
[5] Lu X., Tervola P., Viljanen M., Transient analytical solution to heat conduction in composite circular cylinder, International Journal of Heat and Mass Transfer 2006, 49, 341-348.
[6] Nezhad Y.R., Asemi K., Akhlaghi M., Transient solution of temperature field in functionally graded hollow cylinder with finite length using multi layered approach, International Journal of Mechanics and Materials in Design 2011, 7, 71-82.
[7] Haji-Sheikh A., Beck J.V., Temperature solution in multi-layer bodies, International Journal of Heat and Mass Transfer 2002, 45, 1865-1877.
