BlogFor BusinessTalk to a counsellor
CoursesIIT ProgramsPlacementsLog in
BlogFor Business
Talk to a counsellor
Challenge

Solving the steady and unsteady 2D heat conduction problem

MMohammad Yasir Haneef8 Sept 2021
MATLAB

AIM:

Solve the 2D heat conduction equation by using the point iterative techniques in MATLAB.

OBJECTIVE:

  • To solve the 2D Heat Conduction Equation for Steady and Unsteady state.
  • Using the iterative solver techniques to solve the heat conduction equation.

GIVEN:

The boundary conditions for the problem are as follows:

Top Boundary = 600 K

Bottom Boundary = 900 K

Left Boundary = 400 K

Right Boundary = 800 K

We will implement the following methods for solving implicit equations.

1. Jacobi

2. Gauss-seidel

3. Successive over-relaxation

Our absolute error criteria is 1e-4

THEORY:

In physics and mathematics the heat equation is a partial differential equation that describes how the distribution of some quantity (such as heat) evolves over time in a solid medium, as it spontaneously flows from places where it is higher towards places where it is lower. It is a special case of the diffusion equation.

Steady and Unsteady State:

If the temperature of a body does not vary with time, it is said to be in Steady State. But if there is an abrupt change in its surface temperature, it attains an equilibrium temperature or a steady state after some period. During this period, the temperature varies with time and body is said to be in the unsteady or transient state. This phenomenon is known as Unsteady or transient heat conduction.

Here we use three type of iterative solvers ,

1) Jacobian Iterative Solver :

In jacobi method , each iteration is done using a set of previous values, even if new values are available in the middle of an iteration.Like when you are solving for n+1 iteration and for i th node you will use values of previous(n) iteration even if you have the updated values of the nodes in current iteration.so this method is called method of simultaneous displacements.

The 2D steady state heat conduction equation for this method is

2) Gauss-Seidel iterative solver :

In Gauss-Seidel, as soon as you find a new value in the current iteration you will use it for finding other values in the current iteration.

Like when you are solving for n+1 iteration and for i th node you will use all the available values upto i-1 node of current iteration(n+1) and use i+1 node values and above values of previous iteration(n).SO we get more converged values.so this method is called method of spontaneous displacements.

The 2D steady state heat conduction equation for this method is

3)Successive Over-Relaxation( SOR) method :

It is a variant of the Gauss-Seidel method for solving a linear system of equations, resulting in faster convergence.You can also use jacobi method .Here an extra term called relaxation factor (omega)is multiplied to gauss seidel generated term to scale the convergence rate and solve the equation faster.To calculate new value we take old value and take a parameter omega and take rate of change.

The 2D steady state heat conduction equation for this method is

for scaling jacobian method,where Tj is the result of jacobian method.

and for scaling Gauss-Seidel method,where Tj is the result of Gauss-Seidel method.

when omega = 0 we get new value is equal to old value which gives wrong result.

when omega goes from 0 to 1 we are switching from no update to an update using gauss siedal or jacobi method.

when omega >1 then we are taking the update from gauss siedal method and scaling up the convergence rate. The value we give is called over relaxition factor.

when omega <1 ,The value we give is called under relaxition factor.This is used when the solution given by gauss siedel is over relaxed we give omega value less than 1 to get the converged solution.

MATLAB CODE:

clear all
close all
clc
%intial Conditions
L=1;
nx=20;
ny=nx;
tolerance=1e-4;

error=999;
solver_number= input('enter_the_solver_number: ')

if solver_number==1
steady_jacobian_iteration_method(L,nx,ny,error,tolerance)
end

if solver_number==2
steady_gauss_seidel_iteration_method(L,nx,ny,error,tolerance)
end

if solver_number==3
steady_gauss_SOR_method(L,nx,ny,error,tolerance)
end
  • Solving the steady state 2D heat conduction problem using iterative solvers over a square domain.
  • Here thre different iterative solver used and number is provided for each of them
  • It is stored in a function which is named respectibely.

Function for steady state jacobian iteration method:

function steady_jacobian_iteration_method(L,nx,ny,error,tolerance)

tic
x=linspace(0,L,nx);
y=linspace(0,L,ny);
dx=L/(nx-1);
dy=L/(ny-1);
% TEMPERATURE matrix
T= 300*ones(nx,ny);

T_L = 400;
T_R = 800;
T_T = 600;
T_B = 900;

T(:,1)=T_L;
T(:,end)=T_R;
T(1,:)=T_T;
T(end,:)=T_B;
%edge refining

T(1,1) = (T_L+T_T)/2;
T(1,end) = (T_T+T_R)/2;
T(end,1) = (T_L+T_B)/2;
T(end,end) = (T_B+T_R)/2;
k=(2*(dx^2+dy^2))/(dx^2*dy^2);
Told=T;
j_iteration_number=1;
%convergence loop
while (error>tolerance)

%nodal loop
for i= 2: (nx-1)
for j= 2:(ny-1)

T(i,j)=(1/k)*((Told(i-1,j)+Told(i+1,j))/dx^2 + (Told(i,j-1)+Told(i,j+1))/(dy^2));

end
end
error=max(max(abs(Told-T))); 
%updating old values
Told=T;
j_iteration_number= j_iteration_number+1;
%creating a plot
figure(1)
[C,h]=contourf(x,fliplr(y),T);
clabel(C,h);
colormap(jet)
colorbar;
time_required=toc
title_text = sprintf('solving 2D steady state heat equation using jacobi iterative solver n j_iteration number=%d ; computation time: %0.3f s', j_iteration_number,time_required),
title(title_text);
xlabel('X axis')
ylabel('Y axis')
end
end
  • The function stores the values for the variables given respectively
  • Temperature matrix is stored with respect to nx and ny
  • for loop is used for creating a nodal loop
  • A contour plot is taken fron the values of Temp,r and x.

Function for steady state gauss seidel method:

function steady_gauss_seidel_iteration_method(L,nx,ny,error,tolerance)

tic
x=linspace(0,L,nx);
y=linspace(0,L,ny);
dx=L/(nx-1);
dy=L/(ny-1);
T= 300*ones(nx,ny);

T_L = 400;
T_R = 800;
T_T = 600;
T_B = 900;

T(:,1)=T_L;
T(:,end)=T_R;
T(1,:)=T_T;
T(end,:)=T_B;
%edge refining

T(1,1) = (T_L+T_T)/2;
T(1,end) = (T_T+T_R)/2;
T(end,1) = (T_L+T_B)/2;
T(end,end) = (T_B+T_R)/2;
k=(2*(dx^2+dy^2))/(dx^2*dy^2);
Told=T;

gs_iteration_number=1;
while (error>tolerance)

for i= 2: (nx-1)
for j= 2:(ny-1)

T(i,j)=(1/k)*((T(i-1,j)+Told(i+1,j))/dx^2 + (T(i,j-1)+Told(i,j+1))/(dy^2));

end
end
error=max(max(abs(Told-T)));
Told=T;
gs_iteration_number= gs_iteration_number+1;
figure(1)
[C,h]=contourf(x,fliplr(y),T);
clabel(C,h)
colormap(jet)
colorbar;
time_required=toc;
title_text = sprintf('solving 2D steady state heat equation using gauss seidel iterative solver n iteration number=%d ; computation time: %0.1f s', gs_iteration_number,time_required)
title(title_text);
xlabel('X axis')
ylabel('Y axis')
end
end

Function for Gauss SOR Method:

% Function to solve 2D steady state heat conduction equation using SOR iterative solver.

function steady_gauss_SOR_method(L,nx,ny,error,tolerance)

tic
x=linspace(0,L,nx);
y=linspace(0,L,ny);
dx=L/(nx-1);
dy=L/(ny-1);
T= 300*ones(nx,ny);

T_L = 400;
T_R = 800;
T_T = 600;
T_B = 900;

T(:,1)=T_L;
T(:,end)=T_R;
T(1,:)=T_T;
T(end,:)=T_B;
%edge refining

T(1,1) = (T_L+T_T)/2;
T(1,end) = (T_T+T_R)/2;
T(end,1) = (T_L+T_B)/2;
T(end,end) = (T_B+T_R)/2;
k=(2*(dx^2+dy^2))/(dx^2*dy^2);
Told=T;
omega=1.72;
sor_iteration_number=1;
while (error>tolerance)

for i= 2: (nx-1)
for j= 2:(ny-1)

Tgs(i,j)=(1/k)*((T(i-1,j)+Told(i+1,j))/dx^2 + (T(i,j-1)+Told(i,j+1))/(dy^2));
T(i,j)=(1-omega)*(Told(i,j))+omega*(Tgs(i,j));
end
end
error=max(max(abs(Told-T)));
Told=T;
sor_iteration_number= sor_iteration_number+1;
figure(1)
[C,h]=contourf(x,fliplr(y),T);
clabel(C,h)
colormap(jet)
colorbar;
time_required=toc;
title_text = sprintf('solving 2D steady state heat equation using SOR iterative solver niteration number=%d ; computation time: %0.1f s', sor_iteration_number,time_required);
title(title_text);
xlabel('X axis')
ylabel('Y axis')
end
end

OUTPUTS:

FOR STEADY STATE:

jacobi method

Gauss seidel Method:

SOR Method:

Solving 2D heat conduction equation of transient state :

The equation for 2D heat conduction equation of transient state is

where alpha is the thermal diffusivity.

Then discretize the equation by using forward difference scheme for time and central difference for space.

Here we have a classification,

1) If we consider RHS terms of n th time step while calculating for (n+1) time step then we call it as explicit form.

and simplifying the equation further we get,

As it is of explict form it easy to solve these equations as we use previous time step values which are already known . we solve these equations directly by keeping it in the time loop and space loop.

MATLAB CODE:

%Solving the transient state 2D heat conduction problem using explicit method over a square domain.

clear all 
close all
clc
tic
%intial Conditions
L=1;
nx=20;
ny=nx;
% TEMPERATURE matrix
x=linspace(0,L,nx);
y=linspace(0,L,ny);
dx=L/(nx-1);
dy=L/(ny-1);
% TEMPERATURE matrix
T= 300*ones(nx,ny);

T_L = 400;
T_R = 800;
T_T = 600;
T_B = 900;

T(:,1)=T_L;
T(:,end)=T_R;
T(1,:)=T_T;
T(end,:)=T_B;
%edge refining 

T(1,1) = (T_L+T_T)/2;
T(1,end) = (T_T+T_R)/2;
T(end,1) = (T_L+T_B)/2;
T(end,end) = (T_B+T_R)/2;
Tprev=T;
dt=1e-4;
t=0.15;
nt=t/dt;
alpha=1;
k1=alpha*(dt/dx^2);
k2=alpha*(dt/dy^2);
 
for k=1:nt
    for i=2:nx-1
    for j=2:ny-1
   
    T(i,j)=Tprev(i,j)*(1-2*k1-2*k2)+k1*(Tprev(i-1,j)+Tprev(i+1,j))+k2*(Tprev(i,j-1)+Tprev(i,j+1));
    end
    end
    Tprev=T;
       figure(1)
[C,h]=contourf(x,fliplr(y),T);
clabel(C,h)
colormap(jet)
colorbar;
time_of_simulation=toc;
title_text = sprintf('2D Unsteady state heat equation explicitly n time of simulation=%0.3f',time_of_simulation);
title(title_text)
xlabel('X axis')
ylabel('Y axis')

end
  • The above Code is to solve the 2D Unsteday heat conduction.
  • For loop is compiled to run for different nodal numbers respectively with respect to the given temperature matrix
  • The contour plot is taken to show the exact simulation of the heat conduction with respect to the time of the simulation
  • The title text and labels are given respectively

OUTPUT:

RESULTS:

  • The steady and Unsteady Heat conduction is solved using iterative techniques.
  • The Steady state head conductionn took Approximately 820+ iterations and 400+ iterations with respect to time steps.
  • The unsteady heat conductio took 425 iterations approximately.
  • Among the threee solvers the SOR has taken less time to caluclate with 59 iterations.
Your portfolio, next

Ship work like this
in your own course.

Talk to a career counsellor

Free 20-minute call. We'll map the right course to where you want to work.

By submitting you agree to our Terms and Privacy Policy.

Talk to a career counsellor

Free 20-minute call. We'll map the right course to where you want to work.

By submitting you agree to our Terms and Privacy Policy.