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

Week 5.1 - Mid term project - Solving the steady and unsteady 2D heat conduction problem

RRampally Venkata Sai Kishore16 Jun 2022

Aim : The main aim of the challange is to solve stedy and unstedy state 2D heat conduction equation by using the iterative techniques such as Jacobi method , Gauss-seidal method and Successive over-relaxation method for both implicit and explicit scheme using matlab.

The assumption for solving the method :

  1. Assume that the domain is to be a square.
  2. Assume nx = ny [Number of points along the x directtion is equal to the number of points along y direction].
  3. We are assuming nx = ny as 10.
  4. Boundary conditions for steady and unstedy / transit case, Left = 400 , Right = 800 , Top = 600 , Bottom = 900
  5. Apart from the boundary conditions, all other nodes are given an intial temperature of 300k.
  6. We implement the following methods.
  • Jacobi method
  • Gauss-seidal meythod
  • Sussessive over relaxation method

Theory :

Consider general three dimesion heat comnduction equation which has both time derivative and the spatial terms.

`delT/delt+alpha(del^2T/delx^2 + del^2T/dely^2+del^2T/delz^2)=0`

The above equation has two terms where, term 1 represents time derivative and time 2 is spatial derivative and where alpha is thermal diffusivity.

Now by consedering the assumptions that the above equation is solved for a thin plate, dimensions inthe z dirrection are neglected, thereby the above equation is converts into the following form which is as follows.\

`delT/delt+alpha(del^2T/delx^2 + del^2T/dely^2)=0`

The above equation can be calculated by considering two approaches :

1. Steady state equation :

Under these conditions, the rate of change of temperature remains 0. Therefore in this condition time, the derivative term is equated to zero.

steady state 2D heat conduction equation :

`del^2T/delx^2 + del^2T/dely^2 = 0`

2. Unstedy state equation :

Under these conditions, both terms are considered therefore the equation is follows as :

`delT/delt+alpha(del^2T/delx^2 + del^2T/dely^2)=0`

Above two equations can be solved in 2 ways by using 3 iterative methods such as jacobi, gauss-seidal and SOR.

  1. Implicit method.
  2. Explicit method.

Solving the stedy-state equation :

Consider two dimensional steady state equation

`del^2T/delx^2 + del^2T/dely^2 = 0`

after discreting above equation :

`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)}`

`K = 2.(dx^2+(dy^2)/dx^2.dy^2)`

JACOBI ITERATION[STEADY-STATE] :

CODE :

% JACOBI ITERATION [STEADY-STATE] :
clear all
close all
clc

% Solving the steady state 2D heatconduction (Jocobi method)

% Input conditions 
% Number of nodes
nx = 10;
ny = nx;

x = linspace(0,1,nx);
y = linspace(0,1,ny);

dx = x(2) - x(1);
dy = dx;

% Absolute error criteria > tolerence
error = 9e9;
tolerence = 1e-4;

% input 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;

% Caluclating 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_B)/2;
T(nx,1) = (T_L + T_B)/2;

% Assigning original values of T
T_initial = T;
T_old = T;

% Calculating of 2D stedy heat conduction equation by Jacobi method
iteration_solver = 1;

if iteration_solver == 1
    jacobi_iteration = 1;
    while(error > tolerence)
        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(IMPLICT) = %d' , jacobi_iteration));
pause(0.03);

PLOT :

GAUSS SEIDAL ITERATION [STEADY-STATE] :

CODE :

% GAUSS SEIDAL IRERATION [STEADY-STATE] :
clear all
close all
clc

% Solving the steady state 2D heatconduction (Gauss-seidal method)

% Input conditions 
% Number of nodes
nx = 10;
ny = nx;

x = linspace(0,1,nx);
y = linspace(0,1,ny);

dx = x(2) - x(1);
dy = dx;

% Absolute error criteria > tolerence
error = 9e9;
tolerence = 1e-4;

% input 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;

% Caluclating 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_B)/2;
T(nx,1) = (T_L + T_B)/2;

% Assigning original values of T
T_initial = T;
T_old = T;

% Calculating of 2D stedy heat conduction equation by gauss-seidal method
iteration_solver = 2;

if iteration_solver == 2
    Gauss_seidal_iteration = 1;
    while(error > tolerence)
        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_old(i,j-1) + T_old(i,j+1));
            end
        end
        error = max(max(abs(T_old - T)));
        T_old = T;
        Gauss_seidal_iteration= Gauss_seidal_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 seidal iterations(IMPLICT) = %d' , Gauss_seidal_iteration));
pause(0.03);

PLOT :

SUCCESSIVE OVER RELAXATION ITERATION[STEDY-STADE] :

CODE :

% SUCCESSIVE OVER RELAXATION ITERATION [STEADY-STATE] :
clear all
close all
clc

% Solving the steady state 2D heatconduction (successive_over_relaxation_iteration_method)

% Input conditions 
% Number of nodes
nx = 10;
ny = nx;

x = linspace(0,1,nx);
y = linspace(0,1,ny);

dx = x(2) - x(1);
dy = dx;

% Absolute error criteria > tolerence
error = 9e9;
tolerence = 1e-4;
omega = 1.2;

% input 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;

% Caluclating 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_B)/2;
T(nx,1) = (T_L + T_B)/2;

% Assigning original values of T
T_initial = T;
T_old = T;

% Calculating of 2D stedy heat conduction equation by gauss-seidal method
iteration_solver = 3;

if iteration_solver == 3
    successive_over_relaxation_iteration = 1;
    while(error > tolerence)
        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 iteration(IMPLICT) = %d' , successive_over_relaxation_iteration));
pause(0.03);

PLOT :

JACOBI ITERATION[UNSTEDY-STATE] :

CODE :

%JACOBI ITERATION [UNSTEADY-STATE] :
clear all
close all
clc

% Solving the unsteady state 2D heatconduction (Jocobi method)

% Input conditions 
% Number of nodes
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 > tolerence
error = 9e9;
tolerence = 1e-4;
dt = 1e-3;

% input 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;

% Caluclating 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_B)/2;
T(nx,1) = (T_L + T_B)/2;

% Assigning original values of T
T_initial = T;
T_old = T;

% Calculating of 2D stedy heat conduction equation by Jacobi method
k1 = 1.1*(dt/(dx^2));
k2 = 1.1*(dt/(dx^2));

jacobi_iteration = 1;
for k = 1:nt
    error = 9e9; 
    while(error > tolerence)
        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(1+i,j));
                V = (T_old(i,j+1)) + (T_old(i,j-1));

                T(i,j) = (T_initial(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_initial = 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 unstedy jacobi iterations(EXPLICT) = %d' , jacobi_iteration));
pause(0.03);
end

PLOT :

GAUSS SEIDAL ITERATION [UNSTEDY-STATE] :

CODE :

%GAUSS SEIDAL ITERATION [UNSTEADY-STATE] :
clear all
close all
clc

% Solving the unsteady state 2D heatconduction (Jocobi method)

% Input conditions 
% Number of nodes
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 > tolerence
error = 9e9;
tolerence = 1e-4;
dt = 1e-3;

% input 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;

% Caluclating 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_B)/2;
T(nx,1) = (T_L + T_B)/2;

% Assigning original values of T
T_initial = T;
T_old = T;

% Calculating of 2D stedy heat conduction equation by Jacobi method
k1 = 1.1*(dt/(dx^2));
k2 = 1.1*(dt/(dx^2));

gauss_seidal_iteration = 1;
for k = 1:nt
    error = 9e9; 
    while(error > tolerence)
        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(1+i,j));
                V = (T_old(i,j+1)) + (T_old(i,j-1));

                T(i,j) = (T_initial(i,j)*Term_1)+(H*Term_2)+(V*Term_3);
            end
        end
        error = max(max(abs(T_old - T)));
        T_old = T;
        gauss_seidal_iteration = gauss_seidal_iteration+ 1;
    end


T_initial = 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 unstedy gauss-seidal iterations(IMPLICT) = %d' , gauss_seidal_iteration));
pause(0.03);
end

PLOT :

SUCCESSIVE OVER RELAXATION ITERATION[UN-STEADY] :

CODE :

% SUCESSIVE RELAXATION ITERATION [UNSTEADY-STATE] :
clear all
close all
clc

% Solving the unsteady state 2D heat conduction 

% Input conditions 
% Number of nodes
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 > tolerence
error = 9e9;
tolerence = 1e-4;
dt = 1e-3;
omega = 1.1;


% input 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;

% Caluclating 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_B)/2;
T(nx,1) = (T_L + T_B)/2;

% Assigning original values of T
T_old = T;
T_initial = T;
T_end = T;

% Calculating of 2D stedy heat conduction equation by successive relaxation method
k1 = 1.1*(dt/(dx^2));
k2 = 1.1*(dt/(dx^2));

successive_over_relaxation_iteration = 1;
for k = 1:nt
    error = 9e9; 
    while(error > tolerence)
        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(1+i,j));
                V = (T_old(i,j+1)) + (T_old(i,j-1));

                T(i,j) = (T_initial(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_initial = 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 unstedy successive over relaxation iterations(IMPLICT) = %d' , successive_over_relaxation_iteration));
pause(0.03);
end

PLOT :

EXPLICIT ITERATION :

CODE :

% EXPLICIT ITERATION [STEADY-STATE] :
clear all
close all
clc

% Solving the unsteady state 2D heatconduction (Explicit method)

% Input conditions 
% Number of nodes
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 > tolerence
error = 9e9;
tolerence = 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;

% Caluclating 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_B)/2;
T(nx,1) = (T_L + T_B)/2;

% Assigning original values of T
T_initial = T;
T_old = T;

% Calculating of 2D stedy heat conduction equation by Jacobi method
iteration_solver = 1;
k1 = 1.1*(dt/(dx^2));
k2 = 1.1*(dt/(dx^2));

if iteration_solver == 1
    steady_explicit = 1;
    for k = 1:nt
        error = 1;
    while(error > tolerence)
        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);

PLOT :

Conclusion :

1.In steady state equation SOR methods yeilds the results very quickly when compare to other two methods. Where as jacobi method takes more time to solve the unknows. So the SOR method is the most efficient when compared to other two mwthods to solve the steady state equation.

2.For transit equations , the time required to solve is more when compared to steady state equation. So, therefore it is recommened to solve the 2D heat conduction equation with stady state condition.

In transit state equation, explicit methods yeild the results quickly when compared to implicit method, but it is recomended under the solution when cfl is less than 0.25.

When transit state equation is solved implicitily SOR takes less time when compared to the other two methods where as gauss seidal is more stable when compare to Jacobi.

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.