Numerical Methods in Thermal Engineering

Third Year Mechanical Power - CUFE

Term Project
Due date: Thursday Dec. 21, 2017
Install any version of Matlab on your computer (2015 or later is preferred).
You should have three files to solve this project as given in the following
table,
File Description Status
OneDUnsHeatCond.m Main program Outline provided.
Needs some additions
(inputs, output profiles,
exact solutions)
Coeff1D.m Function to calculate Outline provided.
the active coefficients Needs most of the
work
TDMA.m Tri-Diagonal Matrix Complete
Algorithm solver

The code should solve the energy equation in the following form,
T  2T As q g '''  L2
  Bi    T T   
    2 Acs k

where
T : temperature in C or K
 : normalized time or Fourier number
 : normalized length
 h L
Bi  : Biot number for the ambient air surrounding the domain=
k
A s : the outer surface area
Acs : cross section area
T  : ambient temperature of the surrounding air

q g ''' : volumetric internal heat generation

Using finite difference, derive the computational molecule and rearrange the
equation to have the following form,
AP T i , j 1  AE T i 1 , j 1  AW T i 1 , j 1  B

where,
A P : the active coefficient for the point of interest
A E : the active coefficient east of the point of interest
AW : the active coefficient west of the point of interest
B : the free term containing old time step temperatures as well as any
constants (source terms) form convection and/or internal heat generation.
The given TDMA solver will take care of the equation in that form to solve
the system of equations at the required time step.
This equation is valid for the interior nodes. For the boundary nodes, you
will have to express the dummy nodes in terms of the other adjacent nodes.
You will have to set different values of active coefficients at the boundary
nodes.
You are required to submit a report which includes the following,
A cover page: your full name, section and B.N. are important
Introduction: Include the following
• Discuss what equation you are solving and derive the finite difference
equation.
• Sketch the general computational molecule
• Write down the active coefficients (AP, AE, AW, and B for the interior
nodes (simply use the molecule)
• Derive the active coefficients for the boundary nodes at the LHS and
RHS for type 1, 2 and 3 boundary conditions.
Results and discussion: Include the following,
For each problem you are required to include the following,
• Validate the model (Explicit, Crank Nicolson, and fully implicit) by
comparing the numerical results to the given exact solutions. Use any
value of your choice for  and a value for  which satisfies a
stable “r” value for different methods. For the steady state problems,
you will reach the steady state solution by setting many time steps or a
single very large time step depending on the method you use (implicit,
explicit or CN). Justify your choice of the time step and time steps.
• Plot the numerical solution against the exact solution (at the same
figure)
• Calculate the overall root mean square error for all nodes using this
equation
NP

 T  T numerical 
2
exact
 1
where NP is the number of points
NP
Appendix: Include the following:
• The completed Coeff1D.m file.
• Put the values of all inputs you applied in the main program for each
problem in tables form.
Problem 1
Transient heat conduction in a plane wall with isolated LHS boundary and
convective RHS boundary.
Required: Follow the procedures mentioned in the results and discussion
t
section and calculate the temperature profiles at  0.4535 and
L2
t
 3.2632 . Select a proper time step that satisfies the stability condition for
L2
each scheme.

To  100 C

T  2  0 C
Qs1  0 Bi 2  1


0 1
Exact solution for problem 1
Transient heat conduction in a plane wall with isolated LHS boundary and
convective RHS boundary
t
T T  2  2 2
 C 1 e L cos   
T 0 T  2
where,
T  2 is the ambient temperature at the RHS boundary
T 0 is the initial temperature
C 1 =1.1191
 =0.8603
x

L

T 2
Qs1  0 h2
Bi 2


Problem 2
A very long fin (Ts1 at the LHS boundary and Ts2=T  at the RHS boundary
with convective surrounding)
Required: Follow the procedures mentioned in the results and discussion
section and calculate the steady state temperature profile. Apply a reasonable
initial temperature profile.
Exact solution for problem 2
A very long fin (Ts1 at the LHS boundary and Ts2=T  at the RHS boundary
with convective surrounding)
 Bi  AS 
T T     ACS 

e
T s1 T 
T s1 is the specified temperature at the LHS boundary.
x
 is the normalized length
L
h L
Bi   is the Biot number for the outer surface
k
AS is the outer surface area
ACS is the cross section area
Problem 3
Adiabatic fin (TS1 at the LHS boundary and QS2=0 at the RHS boundary)
and convective surrounding.
Required: Follow the procedures mentioned in the results and discussion
section and calculate the steady state temperature profile. Apply a reasonable
initial temperature profile.
Exact solution for problem 3
Adiabatic fin (TS1 at the LHS boundary and QS2=0 at the RHS boundary)
and convective surrounding


T T  cosh 1    Bi  AS / ACS


T s1 T  
cosh Bi  AS / ACS 
x
 is the normalized length
L
h L
Bi   is the Biot number for the outer surface
k
AS is the outer surface area
ACS is the cross section area
Problem 4
Convective end fin (TS1 at the LHS boundary and Bi2 is given at RHS
boundary = Bi  ) with convective surrounding
Required: Follow the procedures mentioned in the results and discussion
section and calculate the steady state temperature profile. Apply a reasonable
initial temperature profile.
Exact solution for problem 4
Convective end fin (TS1 at the LHS boundary and Bi2 is given at the RHS
boundary = Bi  ) with convective surrounding

cosh  m L 1      sinh  m L 1    
h
T T 

m k 
T s1  T  cosh  m L  
h
sinh  m L 
m k 
Bi   AS / L 
m
L ACS

x
 is the normalized length
L
h L
Bi   is the Biot number for the outer surface
k
AS is the outer surface area
ACS is the cross section area
Problem 5
Transient heat conduction in a plane wall with internal heat generation and
convective boundaries.
Required: Follow the procedures mentioned in the results and discussion
section and calculate the steady state temperature profile. Apply a reasonable
initial temperature profile.

L=0.1m

T 1  60 C T  2  30 C
q g  50 kW / m 3
h1  240 W / m 2 K h2  280 W / m 2 K
k  26 W / mK
hL h2 L
Bi 1  1 Bi 2 
k k


Exact solution for problem 5
Transient heat conduction in a plane wall with internal heat generation and
convective boundaries.

q gL2
T   2  Bi 1 T1  T 1   T1
2k
where
q gL q gL 2 h1
  T 1  Bi 1 T 1  T  2
h2 2k h2
T1 
h1
 Bi 1  1
h2

T 1 T 2
h1 h2
q g
Bi 1 Bi 2