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

AAkshai A28 May 2022

Objective

In this project, here iam going to solve the two dimensional heat conduction for steady and unsteady states.

In case of steady state, here iam solving the equation by using 'three iterative approaches'.

  • Jacobi method
  • Guass siedal method
  • Sucessive over relaxation method.

In case of unsteady state(Transcient state), here iam going to solve the equation in both 'Explicit and Implicit methods'.

Steady state 2D heat conduction equation

Steady state 2D heat conduction equation is given as,

Discretizing the above equation using central differencing scheme we get;

Jacobi method:

Guass-siedal method:

Sucessive over relaxation method:

Unsteady state 2D heat conduction equation

In explicit form:-

Discretizing the above equation using 'forward differencing scheme' for 'time derivative' and 'central differencing scheme' for 'space term'.

In case of explicit problem, CFLx +CFLy <= 0.5 (for a stable solution).

In implicit form:-

Jacobi method:-

Guass-siedal method:-

Sucessive over relaxation method:-

Matlab code

Steady state:-

clear 
close all
clc

% No of nodes
n = 30;
x = linspace(0,1,n);
y = linspace(0,1,n);

% Grid spacing
dx = x(2)-x(1);
dy = dx;

T = 298*ones(n,n);
% Top boundary temp in kelvin
T(1,:) = 600;
% Bottom boundary temp in kelvin
T(end,:) = 900;
% Left boundary temp in kelvin
T(:,1) = 400;
% Right boundary temp in kelvin
T(:,end) = 800;

% Top-left corner temp in kelvin
T(1,1) = 500;
% Top-Right corner temp in kelvin
T(1,end) = 700;
% Bottom-Left corner temp in kelvin
T(end,1) = 650;
% Bottom-Right corner temp in kelvin
T(end,end) = 850;


Told1 = T;
Told2 = T;
Told3 = T;

T1 = T;
T2 = T;
T3 = T;

jacobi_error = 1;
gauss_siedal_error = 1;
SOR_error = 1;
error_req = 1e-4;

jacobi_iter = 0;
gauss_siedal_iter = 0;
SOR_iter = 0;

w = 2/(1+sin(pi*dx));

for method = 1:3
    tic
    if method == 1
        while jacobi_error > error_req
            for i = 2:(n-1)
                for j = 2:(n-1)
                    T1(i,j) = 0.25*(Told1(i+1,j)+Told1(i-1,j)+Told1(i,j+1)+Told1(i,j-1));
                end
            end
            jacobi_iter = jacobi_iter + 1;
            jacobi_error = max(max(abs(Told1-T1)));
            Told1 = T1;
            time_jc = toc;
        end
    end

figure(1)
contourf(x,y,T1,'ShowText','on')
colormap("jet")
colorbar
xlabel('X Axis')
ylabel('Y Axis')
title(['After' num2str(jacobi_iter) 'Jacobi'])
pause(0.3)


tic
if method == 2
    while gauss_siedal_error > error_req
        for i = 2:(n-1)
            for j = 2:(n-1)
                 T2(i,j) = 0.25*(Told2(i+1,j)+T2(i-1,j)+Told2(i,j+1)+T2(i,j-1));
            end
        end

        gauss_siedal_iter = gauss_siedal_iter + 1;
        gauss_siedal_error = max(max(abs(Told2-T2)));
        Told2 = T2;
        time_gs = toc;
    end
end

figure(2)
contourf(x,y,T2,'ShowText','on')
colormap("jet")
colorbar
xlabel('X Axis')
ylabel('Y Axis')
title(['After' num2str(gauss_siedal_iter) 'Gauss Sidel'])
pause(0.3)


tic
if method == 3
    while SOR_error > error_req
        for i = 2:(n-1)
            for j = 2:(n-1)
                T3(i,j) = (1-w)*Told3(i,j) + w*(0.25*(Told3(i+1,j)+T3(i-1,j)+Told3(i,j+1)+T3(i,j-1)));
            end
        end
        SOR_iter = SOR_iter + 1;
        SOR_error = max(max(abs(Told3-T3)));
        Told3 = T3;
        time_sor = toc;
    end
end

figure(3)
contourf(x,y,T2,'ShowText','on')
colormap("jet")
colorbar
xlabel('X Axis')
ylabel('Y Axis')
title(['After' num2str(SOR_iter) 'SOR'])
pause(0.3)

end

Output:-Plots

Simulation time:-

Matlab code:-

Unsteady state:-Explicit:-

clear 
close all
clc

n = 30;
L = 1;
x = linspace(0,L,n);
y = linspace(0,L,n);

dx = x(2)-x(1);
dy = dx;
dt = 2e-4;
alpha = 1.2;

K1 = alpha*dt/(dx^2);
K2 = alpha*dt/(dy^2);
CFL = K1 + K2;

T = 298*ones(n,n);
T(1,:) = 600;
T(end,:) = 900;
T(:,1) = 400;
T(:,end) = 800;

T(1,1) = 500;
T(1,end) = 700;
T(end,1) = 650;
T(end,end) = 850;

Told = T;
T1 = T;

error_req = 1e-4;
error_mag = 1;
iter = 0;

tic
for nt = 1:1000
    for i = 2:(n-1)
        for j = 2:(n-1)

            T1(i,j) = Told(i,j) + K1*(Told(i-1,j)-2*Told(i,j)+Told(i+1,j)) + K2*(Told(i,j-1)-2*Told(i,j)+Told(i,j+1));
        end
    end

    iter = iter+1;
    error_mag = max(max(abs(Told - T1)));
    Told = T1;
    time_req = toc;
    
end

figure(1)
contourf(x,y,T1,'ShowText','on')
colormap("jet")
colorbar
xlabel('X Axis')
ylabel('Y Axis')
title(['Explicit Method at the end of total' num2str(iter) 'Iteration'])

Output:- plots:-

Simulation time:-

Matlab code:-

Unsteady state:-Implicit-Jacobi:- 

clear 
close all
clc

n = 30;
L = 1;
x = linspace(0,L,n);
y = linspace(0,L,n);

dx = x(2)-x(1);
dy = dx;
nt = 1000;
dt = 0.02;
alpha = 1.2e3;

k1 = alpha*dt/(dx^2);
k2 = alpha*dt/(dy^2);
CFL = k1 + k2;

T = 298*ones(n,n);
T(1,:) = 600;
T(end,:) = 900;
T(:,1) = 400;
T(:,end) = 800;

T(1,1) = 500;
T(1,end) = 700;
T(end,1) = 650;
T(end,end) = 850;

jacobi_iter = 0;
jacobi_error = 9e9;
tol = 1e-4;


Told = T;
T_previous = Told;

term1 = (1+2*k1+2*k2)^(-1);
term2 = (k1*term1);
term3 = (k2*term1);



tic
for k = 1:nt
    while jacobi_error > tol
        for i = 2:(n-1)
            for j = 2:(n-1)
                H = Told(i-1,j)+Told(i+1,j);
                V = Told(i,j-1)+Told(i,j+1);
                T(i,j) = (T_previous(i,j)*term1)+(H*term2)+(V*term3);
            end
        end
         jacobi_iter = jacobi_iter + 1;
         jacobi_error = max(max(abs(Told-T)));
         Told = T;
         time_jc = toc;
    end
    T_previous = T;
end

figure(1)
contourf(x,y,T,'ShowText','on')
colormap("jet")
colorbar
xlabel('X Axis')
ylabel('Y Axis')
title(['Implicit Method at the end of total' num2str(jacobi_iter) 'Iteration'])

Output:-plot:-

Simulation time:-

Matlab code:-

Unsteady state:-Implicit-Gauss-seidel:- 

clear 
close all
clc

n = 30;
L = 1;
x = linspace(0,L,n);
y = linspace(0,L,n);

dx = x(2)-x(1);
dy = dx;
nt = 1000;
dt = 0.02;
alpha = 1.2e3;

k1 = alpha*dt/(dx^2);
k2 = alpha*dt/(dy^2);
CFL = k1 + k2;

T = 298*ones(n,n);
T(1,:) = 600;
T(end,:) = 900;
T(:,1) = 400;
T(:,end) = 800;

T(1,1) = 500;
T(1,end) = 700;
T(end,1) = 650;
T(end,end) = 850;

gauss_iter = 0;
gauss_error = 9e9;
tol = 1e-4;


Told = T;
T_previous = Told;

term1 = (1+2*k1+2*k2)^(-1);
term2 = (k1*term1);
term3 = (k2*term1);



tic
for k = 1:nt
    while gauss_error > tol
        for i = 2:(n-1)
            for j = 2:(n-1)
                H = T(i-1,j)+T(i+1,j);
                V = T(i,j-1)+T(i,j+1);
                T(i,j) = (T_previous(i,j)*term1)+(H*term2)+(V*term3);
            end
        end
         gauss_iter = gauss_iter + 1;
         gauss_error = max(max(abs(Told-T)));
         Told = T;
         time_gs = toc;
    end
    T_previous = Told;
end

figure(1)
contourf(x,y,T,'ShowText','on')
colormap("jet")
colorbar
xlabel('X Axis')
ylabel('Y Axis')
title(['Implicit Method at the end of total' num2str(gauss_iter) 'Iteration'])

Output:-plot:-

Simulation time:-

Matlab code:-

Unsteady state:-Implicit-SOR:- 

clear 
close all
clc

n = 30;
L = 1;
x = linspace(0,L,n);
y = linspace(0,L,n);

dx = x(2)-x(1);
dy = dx;
nt = 1000;
dt = 0.02;
alpha = 1.2e3;

k1 = alpha*dt/(dx^2);
k2 = alpha*dt/(dy^2);
CFL = k1 + k2;

T = 298*ones(n,n);
T(1,:) = 600;
T(end,:) = 900;
T(:,1) = 400;
T(:,end) = 800;

T(1,1) = 500;
T(1,end) = 700;
T(end,1) = 650;
T(end,end) = 850;

SOR_iter = 0;
SOR_error = 9e9;
tol = 1e-4;


Told = T;
T_previous = Told;

term1 = (1+2*k1+2*k2)^(-1);
term2 = (k1*term1);
term3 = (k2*term1);

w = 2/(1+sin(pi*dx));
%w = 1.2;



tic
for k = 1:nt
    while SOR_error > tol
        for i = 2:(n-1)
            for j = 2:(n-1)
                H = T(i-1,j)+T(i+1,j);
                V = T(i,j-1)+T(i,j+1);
                T(i,j) = w*((T_previous(i,j)*term1)+(H*term2)+(V*term3)) - (w-1)*T(i,j);
            end
        end
         SOR_iter = SOR_iter + 1;
         SOR_error = max(max(abs(Told-T)));
         Told = T;
         time_SOR = toc;
    end
    T_previous = T;
end

figure(1)
contourf(x,y,T,'ShowText','on')
colormap("jet")
colorbar
xlabel('X Axis')
ylabel('Y Axis')
title(['Implicit Method at the end of total' num2str(SOR_iter) 'Iteration'])

Output:-plot:-

Simulation time:-

Explanation:-

  • Here at first, the number of nodes(n), boundary conditions, x,dx,y,dy and dt values are assigned.
  • The constant such as k for steady state and alpha(k1, k2) ie; the CFL number for transient analysis are defined.
  • For steady state:- Here in case of steady state the solution is analysed by using the ierative solvers such as jacobi,Gauss seidal, SOR. Here the initial errors and the tolerance value is defined and similarly the value of omega is calculated for Sucessive over relaxaion method.
  • Here all iterative techniques are implemented under a for loop. Based on this convergence loop(ie; while loop) is executed.
  • Here Ti,j is calulated under for loop and after the for loop, the error values and iteration values are updated and the value of error is calculated.
  • Tic and Toc commands is used here to calculate the simulation run time.
  • While loop comes to an end when the error falls below the tolerance value and hence Ti,j is plotted using contour and colormap.
  • For transcient state explicit case:-Here in this case time loop for nt = 1:1000 is used and after the loop ends the temperature Ti,j value is calculated and after which the error value and old temperature value is updated. After the Ti,j is plotted at the end of the loop.
  • For transcient state implicit case:- Here in this case the old temperaure value (Told), T_previous and the iteration numbers are assigned. here the convergence loop (ie; while loop) is defined such that Ti,j is calculated under for loop and after the for loop, the error value and iteration values are updated and the error value is calculated.
  • Here Tic and Toc commands is used here to calculate the simulation run time and after every time loop, value of T-previous is updated.
  • While loop comes to an end when the error falls below the tolerance value and hence Ti,j is plotted using contour and colormap.

Here from the above outputs iiam concluding that the jacobi method should take lots of terations to converge, where as the sucessive over relaxation method should takes the least no of iterations to converge as compared that of jacobi and gauss-seidel method. Here the gauss-seidel method lies betwwen the jacobi and sucessive over relaxation method.

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.