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

KKushal Koyyana6 Jul 2023
MATLAB

Objective:

1. Steady-state analysis & Transient State Analysis

Solve the 2D heat conduction equation by using the point iterative techniques that were taught in the class. 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

You will implement the following methods for solving implicit equations.

1. Jacobi

2. Gauss-seidel

3. Successive over-relaxation

Your absolute error criteria are 1e-4.

Solution:

Here, we'll learn about the 2D heat conduction equation's stead-state and transient flow concepts. Before moving on to the coding phase, we must comprehend the derivation phase when the terms are discretized consecutively. Before going over it, there are a few things to keep in mind. For steady-state, we must use iterative solvers, and for transient, we must solve both of the implicit and explicit schemes from the previous part.

Now we are going to understand the derivation part of the steady-state condition.

Solving this by

1. Jacobi

2. Gauss-seidel

3. Successive Over Relaxation SOR

since d2t/dx2 + d2T/dy2 = 0, its has elliptic nature properties which we an discretize them using central differencing scheme

What happens if you pose as Jacob now?

for Gauss seidel

and SOR is an improved Gauss-Seidel method in which the nth and n+1th time steps are identical.

where the relaxation factor, omega, is equal to 1.4.

Here we start the coding portion in MATLAB, where these terms are all taken into account using for loops and while loops.

Program of steady state 2D heat conduction equation.

% Solving Steady and unsteady 2D heat conduction problem
% Steady state conduction 
% Iterative solver - Jacobi, Gauss -seidel, Relaxation methods

clear all
close all
clc

%initial conditions
%L = 1
nx = 10;
ny = 10;
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2) - x(1);
dy = y(2) - y(1);


% given data of boundary condition
% initial temp condition

t = 300*ones(nx,ny);  
t(:,1) = 400;
t(1,:) = 600;
t(:,ny) = 800;
t(nx,:) = 900;

% defining the temperature at corner points

t(1,1) = 500;
t(nx,ny) = 850;
t(1,ny) = 700;
t(nx,1) = 650;

dt = 1e-2; %time step
told=t;

error = 100;
tol = 1e-4;
k = 2*(dx^2+dy^2)/(dx^2*dy^2);
% k2 = (dy^2)/(2*(dx^2+dy^2));

method = input('Iterative method = '); %iterative method selection

%JACOBI ITERATIVE METHOD

if (method==1)
    tic;
    jacobi = 1;
    while(error>tol)
        for i = 2:(nx-1)
            for j = 2:(ny-1)
            t(i,j) = (((told(i-1,j)+told(i+1,j))/(dx^2))+((told(i,j+1)+told(i,j-1))/(dy^2)))/k;
            end
        end
    error = max(max(abs(t-told)));
    told = t;
    jacobi = jacobi+1;
    time_taken = toc;
    end

%plotting    
figure(1)
contourf(x,y,t,'ShowText','on');
colormap(jet);
colorbar;
xlabel('X')
ylabel('Y')
text1 = sprintf('number of iteration = %d', jacobi);
text2 = sprintf('simulation time = %f', time_taken);
title({'Steady state heat condustion profile by jacobi'; text1; text2});
end


%GAUSS-SEIDEL ITERATIVE METHOD

if (method==2)
    tic;
    gauss = 1;
    while(error>tol)
        for i = 2:(nx-1)
            for j = 2:(ny-1)
            t(i,j) = (((t(i-1,j)+told(i+1,j))/(dx^2))+((told(i,j+1)+t(i,j-1))/(dy^2)))/k;
            end
        end
    error = max(max(abs(t-told)));
    told = t;
    gauss = gauss+1;
    time_taken = toc;
    end

%plotting
figure(2)
contourf(x,y,t,'ShowText','on');
colormap(jet);
colorbar;
xlabel('X')
ylabel('Y')
text1 = sprintf('number of iteration = %d', gauss);
text2 = sprintf('simulation time = %f', time_taken);
title({'Steady state heat condustion profile by Gauss-Seidel'; text1; text2});
end


%SUCCESSIVE OVER RELAXATION ITERATIVE METHOD

if (method==3)
    tic;
    sor = 1;
    omega = 1.4; %over relaxation factor
    while(error>tol)
        for i = 2:(nx-1)
            for j = 2:(ny-1)
            term1 =(t(i-1,j)+told(i+1,j))/(dx^2);
            term2 =(told(i,j+1)+t(i,j-1))/(dy^2);
            t(i,j) = ((told(i,j)*(1-omega)))+(omega*(term1+term2))/k;
            end
        end
    error = max(max(abs(t-told)));
    told = t;
    sor = sor+1;
    time_taken = toc;
    end

%plotting
figure(3)
contourf(x,y,t,'ShowText','on');
colormap(jet);
colorbar;
xlabel('X')
ylabel('Y')
text1 = sprintf('number of iteration = %d', sor);
text2 = sprintf('simulation time = %f', time_taken);
text3 = sprintf('over relaxation factor = %f', omega);
title({'Steady state heat condustion profile by Successive Over Relaxation'; text1; text2; text3});
end

Output:

Iterative method By Jacobi Method

Iterative method By Gauss-seidel Method

Iterative method By Successive Over Relaxation method

From the results above, it is clear that the SOR approach performs better in terms of the number of simulations performed and the number of iterations.

Transient (Unsteady) state simulation:

Moving on to the Transient State simulation where we are going to work on Implicit and explicit schemes.

1. Implicit schemes

2. Explicit schemes

Explicit and implicit methods are approaches used in numerical analysis for obtaining numerical approximations to the solutions of time-dependent ordinary and PDEs, as is required in computer simulations of physical process. Explicit methods calculate the state of a system at a later time from the state of the system at the current time, while implicit methods find a solution by solving an equation involving both the current state of the system and the later one. Mathematically,

if

is the current system state and

is the state at the later time

is a small time step), then, for an explicit method

while for an implicit method one solves an equation

By derivation we have here:

Applying the first heat conduction equation in 1 to node at the time moment of , the equation can be rewritten as

The difference quotient can be used to approximate the partial differential in the two sides of (2). The first-order time forward difference quotient can be used to approximate the temperature item to the right of the equal sign.

The second-order partial differential in the left of the equal sign can be approximated by the central difference quotient.

Substituting the value, we can get the difference equation of

The truncation error [34] for equation, which is the difference equation of the heat conduction equation, is O(DT + x2. y2).

Using the most recent equation and the assumption that deltax = deltay = delta, we can derive

While we discretize the whole equation we get here.

By substituting the discrtize values to Jacobi, Gauss-seidel and SOR. we get

Jacobi

Gauss-seidel

Over Relaxation

See transient state below. Implicitly, the original MATLAB code read as follows:

% Solving Steady and Transient 2D heat conduction problem
% Transient state solution with IMPLICIT
% Iterative solver - Jacobi, Gauss -seidel, Relaxation methods

close all
clc

%initial conditions
%L = 1
nx = 10;
ny = 10;
nt = 300; % No.of steps
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2) - x(1);
dy = y(2) - y(1);


% given data of boundary condition
% initial temp condition

t = 300*ones(nx,ny);  
t(:,1) = 400;
t(1,:) = 600;
t(:,ny) = 800;
t(nx,:) = 900;

% defining the temperature at corner points

t(1,1) = 500;
t(nx,ny) = 850;
t(1,ny) = 700;
t(nx,1) = 650;

dt = 0.01; %time step
told=t;
alpha = 0.3; % Thermal diffusivity
tprev_dt = told;

tol = 1e-4;
k1 = alpha*dt/(dx^2);
k2 = alpha*dt/(dy^2);

method = input('Iterative method = '); %iterative method selection

%JACOBI ITERATIVE METHOD

if (method==1)
    tic;
    jacobi = 1;
    for p = 1:nt
        error =100;
         while(error>tol)
        for i = 2:(nx-1)
            for j = 2:(ny-1)
                term1 = (1+(2*k1)+(2*k2))^-1;
                term2 = k1*term1;
                term3 = k2*term1;
                h = (told(i-1,j)+told(i+1,j));
                v = (told(i,j+1)+told(i,j-1));
                t(i,j) = ((tprev_dt(i,j)*term1)+(h*term2)+(v*term3));
            end
        end
        
    error = max(max(abs(t-told)));
    told = t;
    jacobi = jacobi+1;
         end
         tprev_dt = t ;
    end
    time_taken = toc;

%plotting    
figure(1)
contourf(x,y,t,'ShowText','on');
colormap(jet);
colorbar;
xlabel('X')
ylabel('Y')
text1 = sprintf('number of iteration = %d', jacobi);
text2 = sprintf('simulation time = %f', time_taken);
title({'Transient state heat conduction profile Implicit-Jacobi'; text1; text2});
end


%GAUSS-SEIDEL ITERATIVE METHOD

if (method==2)
    tic;
    gauss = 1;
    for p = 1:nt
        error =100;
         while(error>tol)
        for i = 2:(nx-1)
            for j = 2:(ny-1)
                term1 = (1+(2*k1)+(2*k2))^-1;
                term2 = k1*term1;
                term3 = k2*term1;
                h = (t(i-1,j)+told(i+1,j));
                v = (told(i,j+1)+t(i,j-1));
                t(i,j) = ((tprev_dt(i,j)*term1)+(h*term2)+(v*term3));
            end
        end
        
    error = max(max(abs(t-told)));
    told = t;
    gauss = gauss+1;
         end
         tprev_dt = t; 
    end
    time_taken = toc;

%plotting    
figure(1)
contourf(x,y,t,'ShowText','on');
colormap(jet);
colorbar;
xlabel('X')
ylabel('Y')
text1 = sprintf('number of iteration = %d', gauss);
text2 = sprintf('simulation time = %f', time_taken);
title({'Transient state heat conduction profile Implicit-Gauss-seidel'; text1; text2});
end

%SUCCESSIVE OVER RELAXATION ITERATIVE METHOD

if (method==3)
    tic;
    sor = 1;
    omega = 1.2; % Over relaxation factor
    for p = 1:nt
        error =100;
         while(error>tol)
        for i = 2:(nx-1)
            for j = 2:(ny-1)
                term1 = (1+(2*k1)+(2*k2))^-1;
                term2 = k1*term1;
                term3 = k2*term1;
                h = (t(i-1,j)+told(i+1,j));
                v = (told(i,j+1)+t(i,j-1));
                t(i,j) = (told(i,j)*(1-omega))+(omega*((tprev_dt(i,j)*term1)+(h*term2)+(v*term3)));
            end
        end
        
    error = max(max(abs(t-told)));
    told = t;
    sor = sor+1;
         end
         tprev_dt = t ;
    end
    time_taken = toc;

%plotting    
figure(1)
contourf(x,y,t,'ShowText','on');
colormap(jet);
colorbar;
xlabel('X')
ylabel('Y')
text1 = sprintf('number of iteration = %d', sor);
text2 = sprintf('simulation time = %f', time_taken);
title({'Transient state heat conduction profile Implicit-SOR'; text1; text2});
end

Output:

Iterative method = 1

Iterative method = 2

Iterative method = 3

The conclusion of the before mentioned results is identical to that of the steady state, and it has been found that SOR approaches perform better when there are fewer simulations and less run time.

Transient State with Explicit Scheme

after discretizing the values for explicit schemes we have here.

Original Code for Explicit schemes of Transient state 

% Solving Steady and Transient 2D heat conduction problem
% Transient state solution with EXPLICIT

close all
clc

%initial conditions
%L = 1
nx = 10;
ny = 10;
nt = 300; % No.of steps
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2) - x(1);
dy = y(2) - y(1);

% initial temp condition

t = 300*ones(nx,ny);  
t(:,1) = 400;
t(1,:) = 600;
t(:,ny) = 800;
t(nx,:) = 900;

% defining the temperature at corner points

t(1,1) = 500;
t(nx,ny) = 850;
t(1,ny) = 700;
t(nx,1) = 650;

dt = 0.01; %time step
told=t;
alpha = 0.3; % Thermal diffusivity

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

%EXPLICIT METHOD
% running for loop

tic;
    explicit = 1;
    for p = 1:nt
        for i = 2:(nx-1)
            for j = 2:(ny-1)
                term1 = (told(i-1,j)-2*told(i,j)+told(i+1,j));
                term2 = (told(i,j-1)-2*told(i,j)+told(i,j+1));
                t(i,j) = ((told(i,j))+(k1*term1)+(k2*term2));
            end
        end
        
    told = t;
    explicit = explicit+1;
    end
    
  time_taken = toc; 
  
 %Plotting   
figure(1)
contourf(x,y,t,'ShowText','on');
colormap(jet);
colorbar;
xlabel('X')
ylabel('Y')
text1 = sprintf('timen = %d', dt*nt);
text2 = sprintf('simulation time = %f', time_taken);
title({'Transient state heat conduction profile Explicit'; text1; text2});

Output:

Conclusion:

  • The findings of steady state methods converge within a few number of iterations because there is no time loop involved. For implicit methods, we include a time loop that executes repeatedly.
  • Less iterations and computation time are required for the SOR to give us the solution.
  • Explicit approaches have higher accuracy and need less computational time to run.
  • Iterative methods were used to solve the 2D heat conduction equation.

Reference: https://www.hindawi.com/journals/complexity/2019/7432138/

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.