All Courses
All Courses
Courses by Software
Courses by Semester
Courses by Domain
Tool-focused Courses
Machine learning
POPULAR COURSES
Success Stories
OBJECTIVE: To solve Steady and Unsteady state 2D Heat Conduction Equation by using the point iterative techniques such as Jacobi, Gauss Seidel & Successive over-relaxation for both implicit and explicit schemes using MATLAB. Problem statement: Computation area is to be a square. Assume nx = ny [No. of…
Kishoremoorthy SP
updated on 27 Apr 2023
OBJECTIVE:
To solve Steady and Unsteady state 2D Heat Conduction Equation by using the point iterative techniques such as Jacobi, Gauss Seidel & Successive over-relaxation for both implicit and explicit schemes using MATLAB.
Problem statement:
By using Point iterative method, solve the 2D Heat Conduction Equation and plot the Temperature Contour for all the three iterative methods namely Jacobi, Gauss-Siedel & SOR for Steady and Unsteady conditions under IMPLICIT scheme also additionally solving the EXPLICIT scheme for Unsteady Condition.
THEORY:
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.
For a function
of three spatial variables
(see Cartesian coordinate system) and the time variable
, the heat equation is
where,
is a real coefficient called the diffusivity of the medium.
Types of Heat Conduction schemes
This transfer of Heat energy may occur under two conditions
Unsteady state conditions are a precursor to steady state conditions. No system exists initially under steady state conditions. Some time must pass, after heat transfer is initiated, before the system reaches steady state. During that period of transition the system is under unsteady state conditions.
Time Loop Logic:
Boundary Condition Map:
Here,
TL = Left boundary
TT = Top boundary
TP= Middle/Point boundary
TB = Bottom boundary
TR = Right boundary
Governing Equations
In order to solve the 2D Heat Conduction problem, the governing equations for respective conditions with respect to cartesian coordinates should be discretized to enter that in a code which is given under MAIN PROGRAM below.
MAIN PROGRAM:
Below there is a sequence of seven individual programs, starting with IMPLICIT Scheme of Steady & Unsteady-state conditions for three iterative methods Jacobi, Gauss-siedel & SOR followed by EXPLICIT Scheme for Unsteady/Transient Conditions.
IMPLICIT Scheme - Steady-state Condition Program:
Steady-state Equation discretization:
In steady-state conduction, the amount of heat entering any region of an object is equal to the amount of heat coming out (if this were not so, the temperature would be rising or falling, as thermal energy was tapped or trapped in a region).
Jacobi Iteration
close
clear
clc
%Solving the Steady state 2D heat conduction (Jacobi Method)
%Input Parameters
%Number of grid points
nx = 10;
ny = nx;
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2) - x(1);
dy = dx;
%Absolute error criteria > tolerance
error = 9e9;
tolerance = 1e-4;
%Defining Boundary conditions
T_L = 400;
T_T = 600;
T_R = 800;
T_B = 900;
T = 300*ones(nx,ny);
T(2:ny - 1, 1) = T_L;
T(2:ny - 1, nx) = T_R;
T(1, 2:nx - 1) = T_T;
T(ny, 2:nx - 1) = T_B;
%Calculating Average temperature at corners
T(1,1) = (T_T + T_L)/2;
T(nx,ny) = (T_R + T_B)/2;
T(1,ny) = (T_T + T_R)/2;
T(nx,1) = (T_L + T_B)/2;
%Assigning orginal values to T
T_intial = T;
T_old = T;
%Calculation of 2D steady heat conduction EQUATION by JACOBI method
iterative_solver = 1;
if iterative_solver == 1
jacobi_iteration = 1;
while(error > tolerance)
for i = 2:nx - 1
for j = 2:ny - 1
T(i,j) = 0.25*(T_old(i-1,j) + T_old(i+1,j) + T_old(i,j-1) + T_old(i,j+1));
end
end
error = max(max(abs(T_old - T)));
T_old = T;
jacobi_iteration = jacobi_iteration + 1;
end
end
%Plotting
figure(1)
contourf(x,y,T)
clabel(contourf(x,y,T))
colorbar
colormap(jet)
set(gca, 'ydir', 'reverse')
xlabel('X-Axis')
ylabel('Y-Axis')
title(sprintf('No. of Jacobi Iterations(IMPLICIT) = %d', jacobi_iteration));
pause(0.03);
grid on
Gauss-Siedel Iteration
close all
clear all
clc
%Solving the Steady state 2D heat conduction (Gauss-seidel Method)
%Input Parameters
%Number of grid points
nx = 10;
ny = nx;
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2) - x(1);
dy = dx;
%Absolute error criteria > tolerance
error = 9e9;
tolerance = 1e-4;
%Defining Boundary conditions
T_L = 400;
T_T = 600;
T_R = 800;
T_B = 900;
T = 300*ones(nx,ny);
T(2:ny - 1, 1) = T_L;
T(2:ny - 1, nx) = T_R;
T(1, 2:nx - 1) = T_T;
T(ny, 2:nx - 1) = T_B;
%Calculating Average temperature at corners
T(1,1) = (T_T + T_L)/2;
T(nx,ny) = (T_R + T_B)/2;
T(1,ny) = (T_T + T_R)/2;
T(nx,1) = (T_L + T_B)/2;
%Assigning orginal values to T
T_intial = T;
T_old = T;
%Calculation of 2D steady heat conduction EQUATION by Gauss-Seidel method
iterative_solver = 2;
if iterative_solver == 2
Gauss_seidel_iteration = 1;
while(error > tolerance)
for i = 2:nx - 1
for j = 2:ny - 1
T(i,j) = 0.25*(T(i-1,j) + T_old(i+1,j) + T(i,j-1) + T_old(i,j+1));
end
end
error = max(max(abs(T_old - T)));
T_old = T;
Gauss_seidel_iteration = Gauss_seidel_iteration + 1;
end
end
%Plotting
figure(1)
contourf(x,y,T)
clabel(contourf(x,y,T))
colorbar
colormap(jet)
set(gca, 'ydir', 'reverse')
xlabel('X-Axis')
ylabel('Y-Axis')
title(sprintf('No. of Gauss-Seidel Iterations(IMPLICIT) = %d', Gauss_seidel_iteration));
pause(0.03);
grid on
Successive over-relaxation
close all
clear all
clc
%Solving the Steady state 2D heat conduction (Successive over-relaxation Method)
%Input Parameters
%Number of grid points
nx = 10;
ny = nx;
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2) - x(1);
dy = dx;
%Absolute error criteria > tolerance
error = 9e9;
tolerance = 1e-4;
omega = 1.2;
%Defining Boundary conditions
T_L = 400;
T_T = 600;
T_R = 800;
T_B = 900;
T = 300*ones(nx,ny);
T(2:ny - 1, 1) = T_L;
T(2:ny - 1, nx) = T_R;
T(1, 2:nx - 1) = T_T;
T(ny, 2:nx - 1) = T_B;
%Calculating Average temperature at corners
T(1,1) = (T_T + T_L)/2;
T(nx,ny) = (T_R + T_B)/2;
T(1,ny) = (T_T + T_R)/2;
T(nx,1) = (T_L + T_B)/2;
%Assigning orginal values to T
T_intial = T;
T_old = T;
%Calculation of 2D steady heat conduction EQUATION by Successive over-relaxation method
iterative_solver = 3;
if iterative_solver == 3
Successive_Over_Relaxation_iteration = 1;
while(error > tolerance)
for i = 2:nx - 1
for j = 2:ny - 1
T(i,j) = omega*(0.25*(T(i-1,j) + T_old(i+1,j) + T(i,j-1) + T_old(i,j+1))) + (1 - omega)*T_old(i,j);
end
end
error = max(max(abs(T_old - T)));
T_old = T;
Successive_Over_Relaxation_iteration = Successive_Over_Relaxation_iteration + 1;
end
end
%Plotting
figure(1)
contourf(x,y,T)
clabel(contourf(x,y,T))
colorbar
colormap(jet)
set(gca, 'ydir', 'reverse')
xlabel('X-Axis')
ylabel('Y-Axis')
title(sprintf('No. of Successive Over-Relaxation Iterations(IMPLICIT) = %d', Successive_Over_Relaxation_iteration));
pause(0.03);
grid on
Output:
IMPLICIT Scheme Unsteady/Transient-state Condition Program:
Unsteady-state Equation discretization:
Non-steady-state situations appear after an imposed change in temperature at a boundary as well as inside of an object, as a result of a new source or sink of heat suddenly introduced within an object, causing temperatures near the source or sink to change in time.
Jacobi Iteration code
close all
clear all
clc
%Solving the Unsteady state 2D heat conduction by Jacobi Method(IMPLICIT SCHEME)
%Input Parameters
%Number of grid points
nx = 10;
ny = nx;
nt = 1400;
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2) - x(1);
dy = dx;
%Absolute error criteria > tolerance
error = 9e9;
tolerance = 1e-4;
dt = 1e-3;
%Defining Boundary conditions
T_L = 400;
T_T = 600;
T_R = 800;
T_B = 900;
T = 300*ones(nx,ny);
T(2:ny - 1, 1) = T_L;
T(2:ny - 1, nx) = T_R;
T(1, 2:nx - 1) = T_T;
T(ny, 2:nx - 1) = T_B;
%Calculating Average temperature at corners
T(1,1) = (T_T + T_L)/2;
T(nx,ny) = (T_R + T_B)/2;
T(1,ny) = (T_T + T_R)/2;
T(nx,1) = (T_L + T_B)/2;
%Assigning orginal values to T
T_old = T;
T_intial = T;
%Calculation of 2D steady heat conduction EQUATION by Successive over-relaxation method
k1 = 1.1*(dt/(dx^2));
k2 = 1.1*(dt/(dy^2));
Jacobi_iteration = 1;
for k = 1:nt
error = 9e9;
while(error > tolerance)
for i = 2:nx - 1
for j = 2:ny - 1
Term_1 = 1/(1 + (2*k1) + (2*k2));
Term_2 = k1*Term_1;
Term_3 = k2*Term_1;
H = (T_old(i-1,j)) + (T_old(i+1,j));
V = (T_old(i,j+1)) + (T_old(i,j-1));
T(i,j) = (T_intial(i,j)*Term_1) + (H*Term_2) + (V*Term_3);
end
end
error = max(max(abs(T_old - T)));
T_old = T;
Jacobi_iteration = Jacobi_iteration + 1;
end
T_intial = T;
%Plotting
figure(1)
contourf(x,y,T)
clabel(contourf(x,y,T))
colorbar
colormap(jet)
set(gca, 'ydir', 'reverse')
xlabel('X-Axis')
ylabel('Y-Axis')
title(sprintf('No. of Unsteady Jacobi Iterations(EXPLICIT) = %d', Jacobi_iteration));
pause(0.03)
end
Gauss-Siedel Iteration Code
close all
clear all
clc
%Solving the Unsteady state 2D heat conduction by Gauss Seidel Method(IMPLICIT SCHEME)
%Input Parameters
%Number of grid points
nx = 10;
ny = nx;
nt = 1400;
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2) - x(1);
dy = dx;
%Absolute error criteria > tolerance
error = 9e9;
tolerance = 1e-4;
dt = 1e-3;
%Defining Boundary conditions
T_L = 400;
T_T = 600;
T_R = 800;
T_B = 900;
T = 300*ones(nx,ny);
T(2:ny - 1, 1) = T_L;
T(2:ny - 1, nx) = T_R;
T(1, 2:nx - 1) = T_T;
T(ny, 2:nx - 1) = T_B;
%Calculating Average temperature at corners
T(1,1) = (T_T + T_L)/2;
T(nx,ny) = (T_R + T_B)/2;
T(1,ny) = (T_T + T_R)/2;
T(nx,1) = (T_L + T_B)/2;
%Assigning orginal values to T
T_old = T;
T_intial = T;
%Calculation of 2D steady heat conduction EQUATION by Successive over-relaxation method
k1 = 1.1*(dt/(dx^2));
k2 = 1.1*(dt/(dy^2));
Gauss_Seidel_iteration = 1;
for k = 1:nt
error = 9e9;
while(error > tolerance)
for i = 2:nx - 1
for j = 2:ny - 1
Term_1 = 1/(1 + (2*k1) + (2*k2));
Term_2 = k1*Term_1;
Term_3 = k2*Term_1;
H = (T(i-1,j)) + (T_old(i+1,j));
V = (T(i,j+1)) + (T_old(i,j-1));
T(i,j) = (T_intial(i,j)*Term_1) + (H*Term_2) + (V*Term_3);
end
end
error = max(max(abs(T_old - T)));
T_old = T;
Gauss_Seidel_iteration = Gauss_Seidel_iteration + 1;
end
T_intial = T;
%Plotting
figure(1)
contourf(x,y,T)
clabel(contourf(x,y,T))
colorbar
colormap(jet)
set(gca, 'ydir', 'reverse')
xlabel('X-Axis')
ylabel('Y-Axis')
title(sprintf('No. of Unsteady Gauss-Seidel Iterations(IMPLICIT) = %d', Gauss_Seidel_iteration));
pause(0.03)
end
Successive over-relaxation Code:
close all
clear all
clc
%Solving the Unsteady state 2D heat conduction by Successive over-relaxation Method(IMPLICIT SCHEME)
%Input Parameters
%Number of grid points
nx = 10;
ny = nx;
nt = 1400;
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2) - x(1);
dy = dx;
%Absolute error criteria > tolerance
error = 9e9;
tolerance = 1e-4;
dt = 1e-3;
omega = 1.1;
%Defining Boundary conditions
T_L = 400;
T_T = 600;
T_R = 800;
T_B = 900;
T = 300*ones(nx,ny);
T(2:ny - 1, 1) = T_L;
T(2:ny - 1, nx) = T_R;
T(1, 2:nx - 1) = T_T;
T(ny, 2:nx - 1) = T_B;
%Calculating Average temperature at corners
T(1,1) = (T_T + T_L)/2;
T(nx,ny) = (T_R + T_B)/2;
T(1,ny) = (T_T + T_R)/2;
T(nx,1) = (T_L + T_B)/2;
%Assigning orginal values to T
T_old = T;
T_intial = T;
T_end = T;
%Calculation of 2D steady heat conduction EQUATION by Successive over-relaxation method
k1 = 1.1*(dt/(dx^2));
k2 = 1.1*(dt/(dy^2));
Successive_over_relaxation_iteration = 1;
for k = 1:nt
error = 9e9;
while(error > tolerance)
for i = 2:nx - 1
for j = 2:ny - 1
Term_1 = 1/(1 + (2*k1) + (2*k2));
Term_2 = k1*Term_1;
Term_3 = k2*Term_1;
H = (T(i-1,j)) + (T_old(i+1,j));
V = (T(i,j+1)) + (T_old(i,j-1));
T(i,j) = (T_intial(i,j)*Term_1) + (H*Term_2) + (V*Term_3);
T_end(i,j) = ((1-omega)*T_old(i,j)) + (T(i,j)*omega);
end
end
error = max(max(abs(T_old - T)));
T_old = T_end;
Successive_over_relaxation_iteration = Successive_over_relaxation_iteration + 1;
end
T_intial = T_end;
%Plotting
figure(1)
contourf(x,y,T)
clabel(contourf(x,y,T))
colorbar
colormap(jet)
set(gca, 'ydir', 'reverse')
xlabel('X-Axis')
ylabel('Y-Axis')
title(sprintf('No. of Unsteady SOR Iterations(IMPLICIT) = %d', Successive_over_relaxation_iteration));
pause(0.03)
end
Output:
EXPLICIT Scheme - Unsteady/Transient state Condition Program
Unsteady-state Equation discretization(EXPLICIT):
Explicit Scheme code :
close all
clear all
clc
%Solving the Unsteady state 2D heat conduction (EXPLICIT SCHEME)
%Input Parameters
%Number of grid points
nx = 10;
ny = nx;
nt = 1400;
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2) - x(1);
dy = dx;
%Absolute error criteria > tolerance
error = 9e9;
tolerance = 1e-4;
dt = 1e-3;
omega = 1.1;
%Defining Boundary conditions
T_L = 400;
T_T = 600;
T_R = 800;
T_B = 900;
T = 300*ones(nx,ny);
T(2:ny - 1, 1) = T_L;
T(2:ny - 1, nx) = T_R;
T(1, 2:nx - 1) = T_T;
T(ny, 2:nx - 1) = T_B;
%Calculating Average temperature at corners
T(1,1) = (T_T + T_L)/2;
T(nx,ny) = (T_R + T_B)/2;
T(1,ny) = (T_T + T_R)/2;
T(nx,1) = (T_L + T_B)/2;
%Assigning orginal values to T
T_intial = T;
T_old = T;
%Calculation of 2D steady heat conduction EQUATION by Successive over-relaxation method
iterative_solver = 1;
k1 = 1.1*(dt/(dx^2));
k2 = 1.1*(dt/(dy^2));
if iterative_solver == 1
Steady_Explicit = 1;
for k = 1:nt
error = 1;
while(error > tolerance)
for i = 2:nx - 1
for j = 2:ny - 1
Term_1 = T_old(i,j);
Term_2 = k1*(T_old(i+1,j) - 2*T_old(i,j) + (T_old(i-1,j)));
Term_3 = k2*(T_old(i,j+1) - 2*T_old(i,j) + T_old(i,j-1));
T(i,j) = Term_1 + Term_2 + Term_3;
end
end
error = max(max(abs(T_old - T)));
T_old = T;
Steady_Explicit = Steady_Explicit + 1;
end
end
end
%Plotting
figure(1)
contourf(x,y,T)
clabel(contourf(x,y,T))
colorbar
colormap(jet)
set(gca, 'ydir', 'reverse')
xlabel('X-Axis')
ylabel('Y-Axis')
title(sprintf('No. of Unsteady Iterations(EXPLICIT) = %d', Steady_Explicit));
pause(0.03);
Output:
CONCLUSIONS:
From the above outputs of various methods & comparing their No. of iterations, we can conclude that
Leave a comment
Thanks for choosing to leave a comment. Please keep in mind that all the comments are moderated as per our comment policy, and your email will not be published for privacy reasons. Please leave a personal & meaningful conversation.
Other comments...
Week 4 Challenge : CFD Meshing for BMW car
AIM: FOR THE GIVE MODEL, CHECK AND SOLVE ALL GEOMETRICAL ERRORS ON HALF PORTION AND ASSIGN APPROPRITATE PIDS. PERFORMS MESHING WITH GIVEN TARGET LENGTH AND ELEMENT QUALITY CRITERIA. AFTER MESHING THE HALF MODEL,DO SYMMETRY TO THE OTHER SIDE. PRODECURE: INITIALLY, OPEN THE GIVEN BMW MODEL IN ANSA SOFTWARE.…
20 Oct 2023 11:25 AM IST
Week 12 - Validation studies of Symmetry BC vs Wedge BC in OpenFOAM vs Analytical H.P equation
Aim: employing the symmetry boundary condition to simulate an axis-symmetric laminar flow through the pipe's constant cross-section. Using both symmetry and wedge boundary conditions, simulate the aforementioned angles—10, 25, and 45 degrees—and evaluate your results using HP equations. Introduction: Hagen-Poiseuille's…
04 May 2023 03:14 PM IST
Week 11 - Simulation of Flow through a pipe in OpenFoam
Aim: Simulate axisymmetric flow in a pipe through foam. Objective: Verify the hydrodynamic length using the numerical result Verify a fully developed flow rate profile with its analytical profile Verify the maximum velocity and pressure drop for fully developed flow Post-process Shear Stress and verify wall shear stress…
04 May 2023 03:04 PM IST
Week 9 - FVM Literature Review
AIM To describe the need for interpolation schemes and flux limiters in Finite Volume Method (FVM). OBJECTIVE To study and understand What is Finite Volume Method(FVM) Write down the major differences between FDM & FVM Describe the need for interpolation schemes and flux limiters in FVM INTRODUCTION …
03 May 2023 05:47 AM IST
Related Courses
Skill-Lync offers industry relevant advanced engineering courses for engineering students by partnering with industry experts.
© 2025 Skill-Lync Inc. All Rights Reserved.