Lab Project Modification of Transient Heat Conduction Aditya Gopi CH06B001, Ajay Singh Jadun CH06B002, Anshu Gupta CH06B004 BRIEF THEORY When solid body is subject a temperature gradient due to difference in temperature of solid and the fluid around it transient heat conduction occurs. In the presence of the temperature gradient heat flows between the solid and fluid till the time thermal equilibrium is reached. When this process takes place the are two kinds of resistances, internal conductive resistance and the external convective resistance. The ratio of these resistances is called Biot number. Biot number (Bi) determines the validity of a body behaving as a lumped mass. Bi= Heat transfer by convection(Heat transfer by conduction) =hLck
Where
Lc = characteristic length (solid volume/ area) k = thermal conductivity of the solid.
Note: If Biot number (Bi) is << 1 lumped system approximation can be done. In order to find the heat transfer coefficient of the cylinder, we need the use of Heisler charts. These charts have the following characteristics: Plot of T0,0,t- T∞To- T∞ versus Fo= αtLc2 (Fourier Number) with reciprocal of Biot number as a parameter. α is thermal diffusivity = k/ρC and T(0,0,t) is temperature of the centre of the cylinder at time t. Temperature at some location other than at the centre can be determined from the plot of Tr,z,t- T∞To- T∞ versus Fo= αtLc2
with 1Bi as a parameter. (r,z are the
distances measured in the radial direction and z direction respectively from the centre.)
INTRODUCTION We are modifying the current existing Transient Heat Conduction set-up from analogous to digital. We are putting a data collector card which reads temperature and send it as a digital data into a computer for storage. We plot this temperature data, on the Heisler chart obtained 1
after solving 2D transient heat conduction in a cylinder. We then obtain the Biot number, which is used to get the convective heat transfer coefficient. 2-D Transient Heat equation is solved using CFD. The cylinder can be divided into 4 symmetric parts round its axis. It can be seen that all the four parts are symmetric and hence we need only one quadrant for our calculation. The result can be extended to other three.
Figure 1: Division of the cylinder
2
Figure 2: Planer representation of single quadrant Now any point (i,j) in the above quadrant can be represented as:
Figure 3: Representation of any point in the quadrant. 3
Now, centre difference formula at any point can be written as : ∂2T∂z2=Ti-1,j-2Ti,j+T(i+1,j)(dz)2 ∂2T∂z2=Ti,j+1-2Ti,j+T(i,j-1)(dr)2 1r∂T∂r=Ti,j+1-T(i,j-1)2r(dr) 1α∂T∂t=Ti-1,j1dz2+Ti+1,j1dz2+Ti,j-2dz2-2dr2+Ti,j+11dr2+12rdr+T(i,j-1)1dr2-12rdr
Now we can solve this differential equation with the help of boundary conditions, which are listed below. BOUNDARY CONDITIONS Top side node: -k∂T∂r=h(T-Tf) Ti-1,j=Ti+1,j-2hTi,j-Tfdzk Top corner (n,m): -k∂T∂z=h(T-Tf) and -k∂T∂r=h(T-Tf) Ti,j+1=Ti,j-1-2hTi,j-Tfdrk Right side node: -k∂T∂r=h(T-Tf) Ti,j+1=Ti,j-1-2hTi,j-Tfdrk Top side corner (1,1) : -k∂T∂z=h(T-Tf) Ti-1,j=Ti+1,j-2hTi,j-Tfdzk By symmetry T(i,j-1)= T(i,j+1) Bottom corner (1,n) : -k∂T∂r=h(T-Tf) Ti,j+1=Ti,j-1-2hTi,j-Tfdrk By symmetry T(i+1,j)= T(i-1,j)
Left side node : T(i,j-1) = T(i,j+1) 4
Bottom side node: T(i+1,j) = T(i-1,j) Centre node: T(i+1,j) = T(i-1,j) T(i,j+1) = T(i,j-1)
COEFFICIENT CALCULATION Now, we convert the 2-D array into 1-D and for that we need to calculate proper coefficients for each of the temperature parameter. Let us take an example of 3×3 matrix. a11a12a13a21a22a23a31a32a33
Now, if we convert this into 1-D the matrix will be 1×9 with the elements arranged as follows: [a11 a12 a13 a21 a22 a23 a31 a32 a33 ] So, following the same pattern, we can calculate the coefficient of 1-D matrix element corresponding to 2-D element. The calculation is shown below. Central Node: Coeff [n (i-1) +j, n(i-1)+(j-1) ]
1 dr2-12(j-1)(dr)2
Coeff [n (i-1) +j, n(i-1)+j ] Coeff [n (i-1) +j,
-2dz2-2(dr)2
n(i-1)+(j+1)]
1dr2+12(j-1)(dr)2
Coeff [n (i-1) +j, n(i-2)+j ]
1(dz)2
Coeff [n (i-1) +j, n(i)+j ]
1(dz)2
Top side nodes
i=1, j=2 to
(n-1)
Coeff [n (i-1) +j, n (i-1)+j ]
-2(dz)2-2(dr)2-2hk(dz)
Coeff [n (i-1) +j, n (i)+j]
2(dz)2
Coeff [n (i-1) +j, n (i-1)+(j+1)]
1dr2+12(j)(dr)2
Coeff [n (i-1) +j, n (i-1)+(j-1)]
1(dr)2-12(j)(dr)2
Coeff [n (i-1) +j, n (m-1)+(n+1) ] Left side nodes
i=2 to m-1 ,
2hTfk(dz) j=1
Coeff [n(i-1)+j, n(i-2)+j]
1(dz)2
Coeff [n (i-1) +j, n(i)+j]
1(dz)2
Coeff [n (i-1) +j, n(i-1)+j]
-2(dz)2-2(dr)2
Coeff [n (i-1) +j, n(i-1)+(j+1) ]
2(dr)2 5
Bottom side nodes
i=m j=2
(n-1)
Coeff [n (i-1) +j, n (i-2)+j]
2(dz)2
Coeff [n (i-1) +j, n(i-1)+j]
-2(dz)2-2(dr)2
Coeff [n (i-1) +j,
n(i-1)+(j+1)]
1dr2-12(j-1)(dr)2
Coeff [n (i-1) +j,
n(i-1)+(j-1)]
1dr2-12(j-1)(dr)2
Right side i=2
(m-1) , j=n
Coeff [n (i-1) +j, n(i-2)+j]
1(dz)2
Coeff [n (i-1) +j, n(i)+j ]
1(dz)2
Coeff [n (i-1) +j, n(i-1)+j]
-2(dz)2-2(dr)2-2hk1(dr)+12dr(n-1)
Coeff [n (i-1) +j, n(i-1)+(j-1)]
2(dr)2
Coeff [n (i-1) +j, n(m-1)+(n+1)]
2hTfk1dr+12dr(j-1)
Node (1,1) T’(1,1)α=T(2,1)2(dz)2 + T(1,1)-2(dz)2-2(dr)2-2hk(dz) + T(1,2)2(dr)2 +2hTfk(dz)
Now (1,n) (n,. . . . . .) T’(1,1)α= T2,n2dz2 + T1,n-2dz2-2dr2-2hkdz-2hk1dr+12drn-1+ T(1,n-1)2(dz)2+ 2hTfk(dz) + 2hTfk1(dr)+12dr(n-1) Node (m,1)
n(m-1)+1
T’(m,1)α= T(m-1,1)2(dz)2 + T(m,1)-2(dz)2-2(dr)2 + T(m,2)2(dr)2
Node (m,n) T'(m,n)α=Tm-1,n2dz2 + Tm,n-2dz2-2dr2-2hk1dr+12drn-1+ Tm,n-12dz2+ 2hTfk1(dr)+12dr(n-1) Now the temperature can be calculated by passing the array T to ‘ode45’ to solve the transient heat conduction equation.
FIGURES 6
& PLOTS
7
F igu re 6: Hei sle r’s cha rt
8
Fi gu re 7: He isl er’ s ch art
9
Figure 8: Dimensionless temperature profile for heating.
10
Figure 9: Temperature vs. Height (at const. radial distance).
11
Figure 10: Temperature vs. Radial Distance (at constant height)
Figure 11: Dimensionless Temperature vs. Dimensionless Time.
12
Figure 12: Temperature profile for heating.
13
COOLING TEMPERATURE PROFILE
Figure 13: Dimensionless temperature profile for cooling.
14
Figure 14: Temperature vs Height (at constant radial distance).
15
Figure 15: Temperature vs. Radial Distance (at constant height).
Figure 16: Dimensionless Temperature vs. Dimensionless Time.
16
Figure 17: Temperature Profile for cooling.
CODE OUTPUT 17
*****Output for curveplot.m***** ans = Heisler Chart Ready Do you want to plot your data on Heisler chart? (enter 1 for YES) 1 Heating (hit 1) or cooling (hit 2)? (enter zero to ignore) 1 Please enter the filename with extension which contains data for the heating curve heating.txt How many readings for heating? (enter zero to ignore) 33 ans = Convective heat transfer coefficient for heating =1419.26 ans = Biot number =1.215 ans = Lumped system analysis invalid ans = Data plotted. Check Figure Do you want to plot your data on Heisler chart? (enter 1 for YES) 1 Heating (hit 1) or cooling (hit 2)? (enter zero to ignore) 2 Please enter the filename with extension which contains data for the cooling curve cooling.txt How many readings for cooling?(enter zero to ignore) 87 ans = Convective heat transfer coefficient for cooling =29.2913 ans = Biot number =0.025075 ans = Lumped system analysis valid ans = 18
Data plotted. Check Figure Do you want to plot your data on Heisler's chart? (enter 1 for YES) n >> *****Output for Transient.m***** ans = Calculations Completed ans = Centre Temperature (K) =369.5298 ans = Biot number =1.215 ans = Lumped system analysis invalid Do you want to continue? (Hit 1 for yes)1 ans = Calculations Completed ans = Centre Temperature (K) =308.8151 ans = Biot number =0.025074 ans = Lumped system analysis valid Do you want to continue? (Hit 1 for yes)n >>
CONCLUSION The code gives us the Biot number and convective heat transfer coefficient for both heating and cooling. For the node configuration of 12x6 the results are found to be accurate and the execution time too is very small. 19
As found manually the lumped system analysis is valid only for cooling. In case of heating the Biot number was much larger than 0.1. When the temperature profile for heating is studied using ‘Transient.m’, it is observed (by studying the matrix ‘Tf’) that the temperature varies significantly throughout the body at any instant of time. This proves that lumped system analysis is invalid. The code can further be used to study transient heat conduction with other solids and fluids. This lets one have a better understanding of the situation before performing the experiment.
REFERENCES 1. Christopher A Long, (2001), Essential Heat Transfer, Addison Wesley Longman
(Singapore) Pte. Ltd., pp.50-57. 2. M. P. Heisler, (1947), Temperature Charts for Induction and Constant Temperature
Heating, Trans. ASME 69, pp. 227–36. 3. D.J. Higham and N.J. Higham, (2000), MATLAB Guide. SIAM.
4. Retrieved April 10, 2009, from http://www.faculty.virginia.edu/ribando/modules/OneDTransient/ 5. Retrieved April 7, 2009, from http://www.mathworks.com/access/helpdesk/help/toolbox/pde/index.html?/access/hel pdesk/help/toolbox/pde/&http://www.mathworks.com/access/helpdesk/help/helpdesk. html
20