Int. J. Adv. Appl. Math. and Mech.

2(3) (2015) 119 - 125 (ISSN: 2347-2529)

Journal homepage: www.ijaamm.com

International Journal of Advances in Applied Mathematics and Mechanics

Finite element solution to transient asymmetric heat

conduction in multilayer annulus
Research Article

Shubha Verma1, ∗ , V. S. Kulkarni2 , K. C. Deshmukh3

1
Department of Mathematics, Ballarpur Institute of Technology, Ballarpur -442701, Maharashtra, India.
2
PG Department of Mathematics, University of Mumbai, Mumbai-400098, Maharashtra, India.
3
Department of Mathematics, Rashtrasant Tukadoji Maharaj Nagpur University, Nagpur -440010, Maharashtra, India.

Received 04 December 2014; accepted (in revised version) 11 March 2015

Abstract: This study investigate the heat transfer and thermal stresses for the time-dependent asymmetric heat conduction
in a multilayer annulus by using finite element method. The realistic problem of isotropic multilayer annulus
subjected to heat flux with internal heat generation and initial temperature is solved. The results of temperature
and thermal stress have been computed numerically, illustrated graphically and interpreted technically. These
results may be used in many engineering problem like combustion chamber, electric generator, compressor,
hydraulic pumps etc.

MSC: 35G30 • 35K20 • 74S05 • 80A20 • 80M10

Keywords: Heat conduction • Multilayer annulus • Finite element method • Finite difference method • Thermal stresses

c 2015 IJAAMM all rights reserved.

1. Introduction
In present world engineering applications, multilayer components are widely used due to its advantage of com-
bining physical, mechanical, and thermal properties of different materials. Many applications require a detailed
knowledge of transient heat transfer (ie. temperature and heat-flux distribution) within the multilayer component.
Heat conduction in multilayer solids has a many applications in engineering field such as buildings, industrial
furnaces, nuclear reactors, turbines, rockets, space craft, and high-tech devices like instruments etc. Therefore it
has been a topic of continued research interest. Ozisik [1] contains a detailed review of one-dimensional composite
media, mainly focus on orthogonal expansions by using Greens functions and Laplace transform techniques.
Turner et al. [2] developed the basic idea of finite element method to obtain the approximate solution of a compli-
cated problem by replacing it into a simpler problem rather than obtaining the exact solution. In the finite element
method it is often be possible to improve or refine the approximate solution by spending more computational
effort. Dechaumphai et al. [3] used finite element analysis for predicting temperatures and thermal stresses
of heated products. Noon [4] studied the fully discrete formulation of Galerkin partial artificial diffusion finite
element method for solving 2-D coupled Burgers problem using Crank-Nicholson method for the time variable.
The numerical results are obtained by using MATLAB are compared with the exact solution.

Sun et al. [5] presented of a problem of transient heat conduction in a one-dimensional three-layer composite
slab. The Eigen function expansion solution is compared with a finite difference numerical solution. In the
early nineteenth century the manufacturing industry had to face a critical problem of designing an advanced
product with complex geometries, multi-material and different types of boundary conditions. Yang et al. [6] has

∗ Corresponding author.
E-mail address: shubha.katiyar@gmail.com

119
 Finite element solution to transient asymmetric heat conduction in multilayer annulus

established a new thermal stability test for heat conduction in one dimensional multilayer composite solid that
have internal heat generation at a rate proportional to the interior temperature. Dhawan et al. [7] had carried
comparative study of two dimensional transient heat conduction problems with Galerkin finite element method,
Euler modified method, Crank-Nicholson method etc. and determined the temperature decay in aluminum plate.
Jain et al. [8] presented an analytical double-series solution for the time-dependent asymmetric heat conduction
in a multilayer annulus by Eigen values. Singh [9] developed an analytical solution for one-dimensional time
dependent multilayer heat conduction problems by using separation of variable, finite integral transform and
applied this method in Eigen value problems. Singh et al. [10] presented an analytical double-series solution for the
unsteady heat conduction multilayer cylindrical problem by using Separation of variables method and obtained
transient temperature distribution in radial direction. Assume spatially non-uniform and time-independent
volumetric heat sources in each layer. Singh et al. [11] solved the time-dependent heat conduction problem in a
multilayer annulus by using separation of variables method, finite integral transform method and observed the
transverse or radial Eigen values for the solution in polar coordinate system. This study analyzes nuclear fuel
rod subjected to time-dependent boundaries or heat sources. Patil et al. [12] solve multidimensional unsteady
state heat conduction equation with Dirichlet boundary conditions using finite volume numerical technique. The
efficiency of this technique is tested with known analytical solutions and the numerical results obtained. Chen et
al. [13] obtained the solution of time-fractional partial differential equations in a multi-layer annulus by using the
finite integral transform technique and Laplace transform technique.

Thatoi et al. [14] obtained the solution of one dimensional heat flow problem in steady state by using Finite
Difference Method, Finite Volume Method and Finite Element Method and done the comparative analysis with
desired exact solution. A Matlab program was used to find the numerical solution. Kulkarni et al. [15] obtained
the heat transfer and thermal stress analysis of cylinder due to internal heat generation under steady temperature
conditions using integral transform methods. The internal heat generation is taken as cylindrical surface heat
source in annular region of linear length of cylinder and is situated concentrically inside the cylinder.

This paper deals with the realistic problem of the thermal stresses of isotropic multilayer annulus with initial temper-
ature Ti . The finite element formulation has been developed for the solution of governing heat conduction equation
and thermal stress analysis. The Matlab programming is used to evaluate the temperature and thermal stresses in
the multilayer circular annulus. The results of temperature and thermal stress have been computed numerically,
illustrated graphically and interpreted technically.

2. Formulation of the problem

Fig. 1. Schematic diagram of an n-layer annulus.

Consider an n-layer annulus {r0 ≤ r ≤ rn }, shown in Fig. 1. All the layers are assumed to be isotropic in thermal
properties with perfect thermal contact. Let ki and αi be the thermal conductivity and thermal diffusivity of the
i t h layer and independent on temperature. At initial time, each i t h layer is at a specified temperature fi (r, θ ) and

120
 Shubha Verma et al. / Int. J. Adv. Appl. Math. and Mech. 2(3) (2015) 119 - 125

time-independent heat sources g i (r, θ ) are switched on for t > 0. Both the inner (i = 1, r = r0 ) as well as the outer
(i = n , r = rn ) surfaces of the annulus is 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 [[5], [8], [11]].
For such cases, the temperature at the contact interfaces will not be continuous.
The governing heat conduction equation along with the boundary conditions which occupying the space (for ri −1 ≤
r ≤ ri , 0 ≤ θ ≤ 2π, t > 0, where i = 1, 2.....n ) are given below:

1 ∂ 2 Ti (r, θ , t ) 1 ∂ ∂ Ti (r, θ , t ) g i (r, θ ) 1 ∂ Ti (r, θ , t )

+ (r )+ = (1)
r2 ∂ θ2 r ∂r ∂r ki αi ∂t

Boundary and initial conditions are given below:

* Inner surface of the first layer (for 0 ≤ θ ≤ 2π and t > 0),

∂ T1 (r0 , θ , t )
Ai n + Bi n T1 (r0 , θ , t ) = Ci n (θ ) (2)
∂r
* Outer surface of the n t h layer (for 0 ≤ θ ≤ 2π and t > 0),

∂ Tn (rn , θ , t )
Ao u t + Bo u t Tn (rn , θ , t ) = Co u t (θ ) (3)
∂r
* Periodic boundary conditions equation (for ri −1 ≤ r ≤ ri and t > 0, where i = 1, 2.....n )

Ti (r, θ = 0, t ) = Ti (r, θ = 2π, t ) (4)

∂ Ti (r, θ = 0, t ) ∂ Ti (r, θ = 2π, t )

= (5)
∂θ ∂θ

* Interface of the (i − 1)t h and i t h layer equation (for 0 ≤ θ ≤ 2π and t > 0, where i = 1, 2....n )

Ti (ri −1 , θ , t ) = Ti −1 (ri −1 , θ , t ) (6)

∂ Ti (ri −1 , θ , t ) ∂ Ti −1 (ri −1 , θ , t )
ki = ki −1 (7)
∂r ∂r
* Initial condition (for ri −1 ≤ r ≤ ri , 0 ≤ θ ≤ 2π and t = 0, where i = 1, 2....n )

Ti (r, θ , t = 0) = fi (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 Bi n and Bo u t are simultaneously zero is not
considered. In addition, asymmetric boundary conditions can be applied by choosing dependent Ci n and Co u t .
Furthermore, multiple layers with zero inner radius (r0 = 0) can be simulated by assigning zero values to constants
Bi n and Ci n in Eq. (2). Such problem can be solved analytically using finite element method.

3. Solution methodology
3.1. Galerkin finite element approach
Galerkin finite element method is used to convert partial differential equation into algebraic equation [[16], [17],
[18]]. Applying Galerkin method to obtain the residual equations corresponding to equation (1) is

1 ∂ 2 Ti (r, θ , t ) 1 ∂ ∂ Ti (r, θ , t ) g i (r, θ ) 1 ∂ Ti (r, θ , t )

Z Z Z
{ + (r )+ − }N j (r )r d r d θ d z
V
r 2 ∂θ 2 r ∂r ∂r ki αi ∂t

∀ i = 1, 2, ....n , j = 1, 2, ....M

Observing that for the axisymmetric case the integrand is independent of the z co-ordinate equation becomes

1 ∂ 2 Ti (r, θ , t ) 1 ∂ ∂ Ti (r, θ , t ) g i (r, θ ) 1 ∂ Ti (r, θ , t )

Z Z
z { + (r )+ − }N j (r )r d r d θ
A
r 2 ∂θ 2 r ∂r ∂r ki αi ∂t

121
 Finite element solution to transient asymmetric heat conduction in multilayer annulus

Integrating the first terms with respect to θ and second term with respect to r by parts and rearranging the above
equation reduces to

∫ ∫ ∂ N j ∂ Ti 1 ∂ N j ∂ Ti ∂ Ti
 Z Z Z Z 
1
z A rdrdθ + drdθ + Nj rdrdθ
¨∂ r Z∂ Zr AZ
r ∂θ ∂θ
rn §
αi A
2π §
∂t «
1 ∂ Ti 2π ∂ T rn
Z
1
ª ª
=z gi Nj r d r d θ + Nj dr + r Nj dθ
ki A r =r
r ∂ θ θ =0 θ =0
∂ r r =r0
0

Applying the boundary conditions (2) to (7) the above equation becomes

∂ N j ∂ Ti ∂ Ti
Z Z Z Z 
1
z rdrdθ + Nj rdrdθ
¨A ∂ Z
r ∂r
Z
αi
Z 2πA
∂t
Z 2π «
1 rn r0
=z gi Nj r d r d θ + Nj {Co u t (θ ) − Bo u t Tn (rn , θ , t )}d θ − Nj {Ci n (θ ) − Bi n T1 (r0 , θ , t ) d θ
ki A θ =0
Ao u t θ =0
Ai n

Substitute the temperature distribution for two node linear element in the above equation T (r, t ) =
N1 (r )T1 + N2 (r )T2 = [N ]T T and ∂∂ Tt = {T˙ }, equation reduces to
∂ [N ] ∂ [N ]T
R 2π
[N ][N ]T ∂∂ Tt r d r d θ }
r
{T }r d r d θ + α1i z { k1i
RR RR RR
z{ A ∂r ∂r A
= A
g i [N ]r d r d θ + θ =0
[N ] A onu t {Co u t (θ ) −
R 2π r0
Bo u t [N ] {T }}d θ − θ =0 [N ] A i n {Ci n (θ ) − Bi n [N ]T {T }}d θ
T

Finally above equation can be written in matrix form as

[K ]{T } + [C ]{T˙ } = { fQ } + { fg }

3.2. Discretization of multilayer annulus

The multilayer annulus is discretized into 4 numbers of nodes (M) and 3 numbers of linear two nodes elements (NE)
[8], which is shown in the Fig. 2.

Fig. 2. Discretized model of multilayer circular annulus

All the 4 elements are assembled and the final equation is formed as given below

[K ]{T } + [C ]{T˙ } = { fQ } + { fg }

(T (t +δt )−T (t )
Applying the finite difference method as in [4] and substituting {T˙ } = δt the above equation reduces to

(T (t + δt ) − T (t )
[K ]{T } + [C ]{ } = { fQ } + { fg }
δt
this is simplified as under

T (t + δt ) = [C ]−1 { fQ }δt + [C ]−1 { fg }δt − [C ]−1 [K ]{T }δt + {T }

122
 Shubha Verma et al. / Int. J. Adv. Appl. Math. and Mech. 2(3) (2015) 119 - 125

4. Illustrative example and results

A three-layer annulus (for r0 ≤ r ≤ r3 , 0 ≤ θ ≤ 2π) is initially at a uniform zero temperature [8]. For time t > 0, θ -
dependent heat flux given by

q0 θ 2 (π − θ )2 , 0 ≤ θ ≤ π

n
q (r = r3 , θ ) =
0 , π ≤ θ ≤ 2π

is applied at the outer surface (r = r3 ) while the inner surface (r = r0 ) is maintained isothermal at zero tempera-
ture. This leads to the coefficients A i n = 0, Bi n = 1, A o u t = k3 , Bo u t = 0, Ci n (θ ) = 0, and Co u t (θ ) = q n (r3 , θ ) in the
respective boundary condition equations. There is no volumetric heat generation in any of the layers, i.e. g i (r, θ ) = 0.

To demonstrate the results multilayer annulus is divided into three parts (first layer is cast iron, second layer is
aluminum and third layer is of copper). The numerical experimentation has been carried out for the three materials
as iron, aluminum and copper layered composite multilayer annulus with no internal heat generation. Parameter
values used in this problem are r0 = 0, q0 = 1, length (L ) = 0.03m, height (z ) = 0.01m, thermal conductivity k1 =
w w w −6 m 2 2
72.7( m−K 0 ), k 2 = 204.2( m−K 0 ) and k 3 = 386( m−K 0 ), thermal diffusivity α1 = 20.34 ∗ 10 ( s ), α2 = 84.18 ∗ 10−6 ( ms ) and
2
α3 = 112.34 ∗ 10−6 ( ms ), Young modulus of elasticity E1 = 210(G P a ), E2 = 69(G P a ) and E3 = 117(G P a ), coefficient
of linear thermal expansion β1 = 12 ∗ 10−6 ( K10 ), β2 = 22.2 ∗ 10−6 ( K10 ) and β3 = 16.5 ∗ 10−6 ( K10 ).

4.1. Temperature distribution

Using Matlab programming and solving the above equation one obtains the temperatures at all the nodes at different
time. The numerical value of temperature distribution at different time in the multilayer circular annulus due to
point heat source is shown in the Fig. 3.

Fig. 3. Transient temperature variation in the radial direction at different time

4.2. Stress analysis

The Strain relation is

ε = β δt

Stress-strain relation is

σ = E ε = E β δt

The elements thermal stresses along the r-axis at t = 5s e c are obtained using the temperatures at the nodes which
is in the Fir. 4.

123
 Finite element solution to transient asymmetric heat conduction in multilayer annulus

Fig. 4. Thermal Stress variation in the radial direction at different time

5. Concluding remark
In this paper, an analytical solution to the asymmetric transient heat conduction in a layered annulus is presented.
Each layer can have spatially varying and time-independent volumetric heat source. Inhomogeneous boundary
condition of first, second and third kind can be applied in the radial direction. The proposed solution is applicable
to the three layered metallic structure with inner radius r0 = 0 and heat flux q0 = 1. As the thermal diffusivity and
conductivity will increase, the maximum temperature variation is seems and least thermal stresses in a material is
generated. The results presented in this paper should provide important information and useful guidance for mul-
tilayer annulus in many engineering problem like combustion chamber, electric generator, compressor, hydraulic
pumps etc. in automobiles field.

References

[1] M. N. Ozisik, Heat Conduction, John Wiley and Sons, New York, 1993.
[2] M. J. Turner, R. W. Clough, H. C. Martin and L. J. Topp, Stiffness and deflection analysis of complex structures,
Journal of Aeronautical Sciences 23 (1956) 805-824.
[3] Pramote Dechaumphai and Wiroj Lim, Finite Element Thermal-Structural Analysis of Heated products, Chu-
lalongkorn University Press, Bangkok, 1996.
[4] Najat Jaleel Noon, Fully discrete formulation of Galerkin-Partial artificial diffusion finite element method for
coupled Burgers problem, International Journal of Advances in Applied Mathematics and Mechanics 1(3) (2014)
56-75.
[5] Yuzhi Sun, Indrek S. Wichman, On Transient heat conduction in a one-dimensional composite slab, Interna-
tional Journal of Heat and Mass Transfer 47 (2004) 1555-1559.
[6] Bingen Yang, Hang Shi, A Thermal Stability criterion for Heat conduction in Multilayer Composite Solids, Jour-
nal of Heat Transfer 131 (2009) 1-6.
[7] Sharanjeet Dhawan, Sheo Kumar, A comparative study of numerical techniques for 2D transient heat conduc-
tion equation using Finite Element method, International Journal of Research and Reviews in Applied Science
1 (2009) 38-46.
[8] Prashant K. Jain, Suneet Singh and Rizwan Uddin , Analytical Solution to Transient Asymmetric Heat Conduc-
tion in a Multilayer Annulus, Journal of Heat Transfer 131 (2009) 011304-1-7.
[9] Suneet Singh , Analytical Solution of Time-Dependent Multilayer Heat Conduction Problems for Nuclear Appli-
cations, Proceedings of the 1st International Nuclear and Renewable Energy Conference (INREC10), Amman,
Jordan, March 21-24 2010.
[10] Suneet Singh , Prashant K. Jain and Rizwan Uddin, Analytical solution to transient heat conduction in polar
coordinates with multiple layers in radial direction, International journal of thermal Sciences 47(3) (2008) 261-
273.
[11] Singh Suneet, Jain K. Prashant and Uddin Rizwan, Finite integral transform method to solve asymmetric heat
conduction in a multilayer annulus with time-dependent boundary conditions, Nuclear Engineering and De-
sign 241(1) (2011) 144-154.

124
 Shubha Verma et al. / Int. J. Adv. Appl. Math. and Mech. 2(3) (2015) 119 - 125

[12] Parag V. Patil, J. S. V. R. Krishna Prasad, The unsteady state finite volume numerical grid technique for multi-
dimensional problems, International Journal of Advances in Applied Mathematics and Mechanics 2(2) (2014)
78-87.
[13] S. Chen, X. Jiang, Analytical solutions to time-fractional partial differential equations in a two-dimensional
multilayer annulus, Statistical Mechanics and its Application, 391(15)(2012) 3865-3874.
[14] N. D. Thatoi, A. Acharya, N. Anand, P. Acharya and V. Agrawal, Numerical Solution of Radial Heat Conduction in
an Infinitely Long Cylinder, Proceedings of the International Conference on Mechanical Engineering, Dhaka,
Bangladesh (ICME2011) 18-20 December, 2011.
[15] V. S. Kulkarni, K. C. Deshmukh and P. H. Munjankar, Thermal stress analysis due to surface heat source, Inter-
national Journal of Advances in Applied Mathematics and Mechanics 2(2) (2014) 139-149.
[16] Robert D. Cook , David S. Malkus and Michael E. Plesha, Concept and Application of Finite Element Analysis,
John Wiley and Sons, New York, third edition, 2000.
[17] David Hutton, Fundamentals of Finite Element Analysis, Tata McGraw-Hill Publishing Company Limited, New
York, second edition, 2006.
[18] Singiresu S. Rao, The Finite Element Method in Engineering, Elsevier Publication, U. K., fourth edition, 2005.

125