Solving the steady and transient 2D heat conduction equations.

Objective

Solving the 2D heat conduction equation using different types of iterative solvers. The parameters used are as follows

The flow of information in this report will be as follows

1. Understanding of the problem

2. Main Code

3. Output plots

4. Inference

1. Understanding of the problem

We know the first order 2D heat conduction equation to be

`(∂T)/(∂t) =α[(∂^2 T)/(∂x^2 )+(∂^2 T)/(∂y^2 )]` (1)

Now eq. (1) can be solved for both steady and transient state as follows

(a) Steady State

At steady state, the 2D heat conduction equation is

`[(∂^2 T)/(∂x^2 )+(∂^2 T)/(∂y^2 )]=0` (2)

Eq. (2) is an elliptical PDE, therefore, discretizing eq.(2) using central difference scheme, we get

or, `[T_(i-1,j) - 2T_(i,j) + T_(i+1,j))/(Δx^2 )+(T_(i,j-1) - 2T_(i,j) + T_(i,j+1))/(Δy^2 ) =0`

or, `(T_(i-1,j ) + T_(i+1,j))/(Δx^2 )+(T_(i,j-1 ) + T_(i,j+1))/(Δy^2 ) = (2T_(i,j))/(Δx^2 )+(2T_(i,j))/(Δy^2 )`

or, `(T_(i-1,j ) + T_(i+1,j))/(Δx^2 )+(T_(i,j-1 ) + T_(i,j+1))/(Δy^2 ) = [2/(Δx^2 )+2/(Δy^2 )]T_(i,j)`

Assuming, `[2/(Δx^2 )+2/(Δy^2 )] = k`

or, `(T_(i-1,j ) + T_(i+1,j))/(Δx^2 )+(T_(i,j-1 ) + T_(i,j+1))/(Δy^2 ) = kT_(i,j)`

or, `(1/k)[(T_(i-1,j ) + T_(i+1,j))/(Δx^2 )+(T_(i,j-1 ) + T_(i,j+1))/(Δy^2 )] = T_(i,j)`

Therefore, the discretized form of the steady state 2D heat conduction equation is

`T_(i,j) = (1/k)[(T_(i-1,j ) + T_(i+1,j))/(Δx^2 )+(T_(i,j-1 ) + T_(i,j+1))/(Δy^2 )]` (3)

Now eq.(3) can be solved using three different types of iterative solvers, namely

(i) Jacobi Iterative Solver

(ii) Gauss-Seidel Solver

(iii) Successive-Over Relaxation Solver

Eq.(3) for Jacobi iterative solver becomes

`Tn_(i,j) = (1/k)[(T_(i-1,j ) + T_(i+1,j))/(Δx^2 )+(T_(i,j-1 ) + T_(i,j+1))/(Δy^2 )]`

Eq.(3) for Gauss-Seidel iterative solver becomes

`Tn_(i,j) = (1/k)[(Tn_(i-1,j ) + T_(i+1,j))/(Δx^2 )+(Tn_(i,j-1 ) + T_(i,j+1))/(Δy^2 )]`

For SOR, we will be implementing it on Gauss-Seidel solver, therefore, eq.(4) for SOR on Gauss-Seidel solver becomes

`Tn_(i,j )={(1-ω)T_((i,j)}}+{ω(1/k [(Tn_(i-1,j ) + T_(i+1,j))/(Δx^2 )+(Tn_(i,j-1 ) + T_(i,j+1))/(Δy^2 )])}`

(b) Transient State

Eq. (1) is a parabolic PDE, therefore, the time marching can be discretized using the two types of time integration schemes and the spatial marching can be discretized using central difference scheme,

Explicit time integration for time marching and central difference scheme for spatial marching

`(∂T)/(∂t) =α[(∂^2 T)/(∂x^2 )+(∂^2 T)/(∂y^2 )]`

or, `[(T_(i,j)^(t+1)-T_(i,j)^t)/(Δt)]=α[(T_(i-1,j) - 2T_(i,j) + T_(i+1,j))/(Δx^2 )+(T_(i,j-1) - 2T_(i,j) + T_(i,j+1))/(Δy^2 )]^t`

or, `T_(i,j)^(t+1)=T_(i,j)^t+αΔt[(T_(i-1,j) - 2T_(i,j) + T_(i+1,j))/(Δx^2 )+(T_(i,j-1) - 2T_(i,j) + T_(i,j+1))/(Δy^2 )]^t`

Assuming, `[(αΔt)/(Δx^2 )]=k_1`, and `[(αΔt)/(Δy^2 )]=k_2`

or, `T_(i,j)^(t+1)=T_(i,j)^t+k_1 [T_(i-1,j) - 2T_(i,j) + T_(i+1,j) ]^t+k_2 [T_(i,j-1) - 2T_(i,j) + T_(i,j+1) ]^t`

Further simplifying the above equation gives,

`T_(i,j)^(t+1)=[1 - 2k_1- 2k_2 ]T_(i,j)^t+k_1 [T_(i-1,j)+ T_(i+1,j) ]^t+k_2 [T_(i,j-1) + T_(i,j+1) ]^t` (4)

Now eq.(4) can be solved using three different types of iterative solvers, namely

(i) Jacobi Iterative Solver

(ii) Gauss-Seidel Solver

(iii) Successive-Over Relaxation Solver

Eq.(4) for Jacobi iterative solver becomes

`T_(i,j)^(t+1)=[1 - 2k_1- 2k_2 ]T_(i,j)^t+k_1 [T_(i-1,j)+ T_(i+1,j) ]^t+k_2 [T_(i,j-1) + T_(i,j+1) ]^t`

Eq.(4) for Gauss-Seidel iterative solver becomes

`T_(i,j)^(t+1)=[1 - 2k_1- 2k_2 ]T_(i,j)^t+k_1 [T_(i-1,j)^(t+1)+ T_(i+1,j)^t]+k_2 [T_(i,j-1)^(t+1) + T_(i,j+1)^t]`

For SOR, we will be implementing it on Gauss-Seidel solver, therefore, eq.(4) for SOR on Gauss-Seidel solver becomes

`T_(i,j )={(1-ω)T_((i,j)}}+{ω([1 - 2k_1- 2k_2 ]T_(i,j)^t+k_1 [T_(i-1,j)^(t+1)+ T_(i+1,j)^t]+k_2 [T_(i,j-1)^(t+1) + T_(i,j+1)^t])}`

Implicit time integration for time marching and central difference scheme for spatial marching

`(∂T)/(∂t) =α[(∂^2 T)/(∂x^2 )+(∂^2 T)/(∂y^2 )]`

or, `[(T_(i,j)^(t+1)-T_(i,j)^t)/(Δt)]=α[(T_(i-1,j) - 2T_(i,j) + T_(i+1,j))/(Δx^2 )+(T_(i,j-1) - 2T_(i,j) + T_(i,j+1))/(Δy^2 )]^(t+1)`

or, `T_(i,j)^(t+1)=T_(i,j)^t+αΔt[(T_(i-1,j) - 2T_(i,j) + T_(i+1,j))/(Δx^2 )+(T_(i,j-1) - 2T_(i,j) + T_(i,j+1))/(Δy^2 )]^(t+1)`

Assuming, `[(αΔt)/(Δx^2 )]=k_1`, and `[(αΔt)/(Δy^2 )]=k_2`

or, `T_(i,j)^(t+1)=T_(i,j)^t+k_1 [T_(i-1,j) - 2T_(i,j) + T_(i+1,j) ]^(t+1)+k_2 [T_(i,j-1) - 2T_(i,j) + T_(i,j+1) ]^(t+1)`

Further simplifying the above equation gives,

`T_(i,j)^(t+1)=(T_(i,j)^t+k_1 [T_(i-1,j)+ T_(i+1,j) ]^(t+1)+k_2 [T_(i,j-1) + T_(i,j+1) ]^(t+1))/(1 + 2k_1 + 2k_2)` (5)

Now eq.(4) can be solved using three different types of iterative solvers, namely

(i) Jacobi Iterative Solver

(ii) Gauss-Seidel Solver

(iii) Successive-Over Relaxation Solver

Eq.(5) for Jacobi iterative solver becomes

`Tn_(i,j)^(t+1)=(T_(i,j)^t+k_1 [T_(i-1,j)+ T_(i+1,j) ]^(t+1)+k_2 [T_(i,j-1) + T_(i,j+1) ]^(t+1))/(1 + 2k_1 + 2k_2)`

Eq.(5) for Gauss-Seidel iterative solver becomes

`Tn_(i,j)^(t+1)=(T_(i,j)^t+k_1 [Tn_(i-1,j)+ T_(i+1,j) ]^(t+1)+k_2 [Tn_(i,j-1) + T_(i,j+1) ]^(t+1))/(1 + 2k_1 + 2k_2)`

For SOR, we will be implementing it on Gauss-Seidel solver, therefore, eq.(5) for SOR on Gauss-Seidel solver becomes

`Tn_(i,j )={(1-ω)T_((i,j)}}+{ω((T_(i,j)^t+k_1 [Tn_(i-1,j)+ T_(i+1,j) ]^(t+1)+k_2 [Tn_(i,j-1) + T_(i,j+1) ]^(t+1))/(1 + 2k_1 + 2k_2))}`

2. Main Program

clear all
close all
clc

%% Parameters

L = 1;                          %Length
n = 10;                         %Number of nodes
dx = L/(n-1);                   %Grid size in x-direction
dy = dx;                        %Grid size in y-direction
x = 0:dx:L;                     %Grid creation in x-axis
y = 0:dy:L;                     %Grid creation in y-axis
tolerance_error = 1e-4;         %Absolute error tolerance
new_error = 1e9;                %Absolute guess error
alpha = 0.3;                    %Diffusivity coefficient
omega = 1.5;                    %Relaxation factor
[xx,yy] = meshgrid(x,y);        %Meshgrid for contour plot

%% Case Setup

% Temperature grid
T = ones(n,n);

% Boundary Conditions (Units in K)
T(:,1) = 400;
T(:,end) = 800;
T(1,:) = 600;
T(end,:) = 900;

% Edge refinement
T(1,1) = (400+900)/2;
T(end,1) = T(1,1);
T(end,end) = (800+900)/2;
T(1,end) = T(end,end);

% Temperature grid storage
T_prev_timestep = T;
T_in = T;
T_old = T_in;

%% User Input
solver_type = input('Enter the type of state you want to perform the heat conduction:\n 1. Steady State\n 2. Transient State\n');

%% Steady State
if(solver_type == 1)

% Task 
   iter_type_ss = input('\nEnter the type of iteration you want to perfom:\n 1. Jacobi\n 2. Gauss Seidel\n 3. SOR on Gauss-Siedel\n');
   if (iter_type_ss == 1)
       jacobi_type = jacobi_iter_ss(new_error, tolerance_error, n, xx, yy, dx, dy, T, T_old);
   elseif (iter_type_ss == 2)
       gs_type = gs_iter_ss(new_error, tolerance_error, n, xx, yy, dx, dy, T_in, T_old);
   else
       sor_gs_type = sor_gs_iter_ss(new_error, tolerance_error, n, xx, yy, omega, dx, dy, T_in, T_old);
   end
end

%% Transient State
if(solver_type == 2)

    inte_type = input('\nEnter the type of integration scheme you want to perfom on the transient heat conduction equation:\n 1. Explicit\n 2. Implicit\n');
    if (inte_type == 1)
% Explicit Task 
   iter_type_ts_ex = input('\nEnter the type of iteration you want to perfom:\n 1. Jacobi\n 2. Gauss Seidel\n 3. SOR on Gauss-Siedel\n');
       if (iter_type_ts_ex == 1)
           jacobi_type = jacobi_iter_tse(n, new_error, tolerance_error, alpha, dx, dy, xx, yy, T_in, T_old);
       elseif (iter_type_ts_ex == 2)
           gs_type = gs_iter_tse(n, new_error, tolerance_error, alpha, dx, dy, xx, yy, T_in, T_old);
       else
           sor_gs_type = sor_gs_iter_tse(n, new_error, tolerance_error, alpha, dx, dy, xx, yy, omega, T_in, T_old);
       end
       
% Implicit Task 
      else
        iter_type_ts_im = input('\nEnter the type of iteration you want to perfom:\n 1. Jacobi\n 2. Gauss Seidel\n 3. SOR on Gauss-Siedel\n');
           if (iter_type_ts_im == 1)
               jacobi_type = jacobi_iter_tsi(new_error, tolerance_error, alpha, n, xx, yy,  dx, dy, T_in, T_old, T_prev_timestep);
           elseif (iter_type_ts_im == 2)
               gs_type = gs_iter_tsi(new_error, tolerance_error, alpha, n, xx, yy,  dx, dy, T_in, T_old, T_prev_timestep);
           else
               sor_gs_type = sor_gs_iter_tsi(new_error, tolerance_error, alpha, n, xx, yy,  dx, dy, omega, T_in, T_old, T_prev_timestep);
           end
    end 
end



%% Steady State Functions

% Steady_Jacobi
function T_in = jacobi_iter_ss(new_error, tolerance_error, n, xx, yy, dx, dy, T_in, T_old)
k = ((2/(dx^2))+(2/(dy^2)));
j_iter = 1;
tic;
while(new_error > tolerance_error)
for i = 2:n-1
    for j = 2:n-1
        T_in(i,j) = (((T_old(i-1,j)+T_old(i+1,j))/dx^2) + ((T_old(i,j-1)+T_old(i,j+1))/dy^2))/k;
    end
end
new_error = max(max(abs(T_old-T_in)));
T_old = T_in;
j_iter = j_iter + 1;
end
tj = toc;
[c,d]=contourf(xx,yy,T_in);
colormap(jet)
colorbar
xlabel('X-axis')
ylabel('Y-axis')
clabel(c,d);
colorbar
title(sprintf('Steady state using Jacobi iterative solver\n Number of nodes = %d\n Number of iteration= %d ; Computation Time = %fs\n Error = %f\n', n,j_iter,tj,new_error))
saveas(gca,'Steady_Jacobi.jpeg')
end

% Steady_Gauss Seidel
function T_in = gs_iter_ss(new_error, tolerance_error, n, xx, yy, dx, dy, T_in, T_old)
k = ((2/(dx^2))+(2/(dy^2)));
g_iter = 1;
tic;
while(new_error > tolerance_error)
for i = 2:n-1
    for j = 2:n-1
        T_in(i,j) = (((T_in(i-1,j)+T_old(i+1,j))/dx^2) + ((T_in(i,j-1)+T_old(i,j+1))/dy^2))/k;
    end
end
new_error = max(max(abs(T_old-T_in)));
T_old = T_in;
g_iter = g_iter + 1;    
end
tg = toc;

[c,d] = contourf(xx,yy,T_in);
colormap(jet)
colorbar
xlabel('X-axis')
ylabel('Y-axis')
clabel(c,d);
colorbar
title(sprintf('Steady state using Gauss-Seidel iterative solver\n Number of nodes = %d\n Number of iteration= %d ; Computation Time = %fs\n Error = %f\n', n,g_iter,tg,new_error))
saveas(gca,'Steady_GS.jpeg')
end

% Steady_SOR on Gauss-Seidel
function T_in = sor_gs_iter_ss(new_error, tolerance_error, n, xx, yy, omega, dx, dy, T_in, T_old)
k = ((2/(dx^2))+(2/(dy^2)));
sg_iter = 1;
tic;
while(new_error > tolerance_error)
for i = 2:n-1
    for j = 2:n-1
        T_in(i,j) = ((1-omega)*T_old(i,j)) + (omega*(((T_in(i-1,j)+T_old(i+1,j))/dx^2) + ((T_in(i,j-1)+T_old(i,j+1))/dy^2))/k);
    end
end
new_error = max(max(abs(T_old-T_in)));
T_old = T_in;
sg_iter = sg_iter + 1;  
end
tsg = toc;

[c,d] = contourf(xx,yy,T_in);
colormap(jet)
colorbar
xlabel('X-axis')
ylabel('Y-axis')
clabel(c,d);
colorbar
title(sprintf('Steady state using SOR on Gauss-Seidel iterative solver\n Number of nodes = %d\n Number of iteration= %d ; Computation Time = %fs\n Error = %f\n', n,sg_iter,tsg,new_error))
saveas(gca,'Steady_SOR_GS.jpeg')
end

%% Transient State Functions

% Explicit Functions

% Transient_Explicit_Jacobi
function T_in = jacobi_iter_tse(n, new_error, tolerance_error, alpha, dx, dy, xx, yy, T_in, T_old)
dt = 0.01;
k1 = (alpha * dt)/dx^2;
k2 = (alpha * dt)/dy^2;
tje_iter = 1;
nt = 200;

tic;
for tl = 1:nt
    while(new_error > tolerance_error)
        for i = 2:n-1
            for j = 2:n-1
                 T_in(i,j) = ((1-(2*k1)-(2*k2))*T_old(i,j)) + ((T_old(i-1,j)+T_old(i+1,j))*k1) + ((T_old(i,j-1)+T_old(i,j+1))*k2); 
            end
        end
    new_error = max(max(abs(T_old-T_in)));
    T_old = T_in;
    tje_iter = tje_iter + 1;
    end
end
ttji = toc;

[c,d] = contourf(xx,yy,T_in);
colormap(jet)
colorbar
xlabel('X-axis')
ylabel('Y-axis')
clabel(c,d);
title(sprintf('Transient state (Explicit) using Jacobi iterative solver\n Number of nodes = %d  ; Timestep = %.2f\n Number of iteration= %d ; Computation Time = %fs\n Error = %f\n', n,dt,tje_iter,ttji,new_error))
saveas(gca,'Transient_Jac_ex.jpeg')
end

% Transient_Explicit_Gauss Seidel
function T_in = gs_iter_tse(n, new_error, tolerance_error, alpha, dx, dy, xx, yy, T_in, T_old)
dt = 0.01;
k1 = (alpha * dt)/dx^2;
k2 = (alpha * dt)/dy^2;
tge_iter = 1;
nt = 200;

tic;
for tl = 1:nt
    while(new_error > tolerance_error)
        for i = 2:n-1
            for j = 2:n-1
                 T_in(i,j) = ((1-(2*k1)-(2*k2))*T_old(i,j)) + ((T_in(i-1,j)+T_old(i+1,j))*k1) + ((T_in(i,j-1)+T_old(i,j+1))*k2); 
            end
        end
    new_error = max(max(abs(T_old-T_in)));
    T_old = T_in;
    tge_iter = tge_iter + 1;
    end
end
ttge = toc;

[c,d] = contourf(xx,yy,T_in);
colormap(jet)
colorbar
xlabel('X-axis')
ylabel('Y-axis')
clabel(c,d);
colorbar
title(sprintf('Transient state (Explicit) using Gauss Seidel iterative solver\n Number of nodes = %d  ; Timestep = %.2f\n Number of iteration= %d ; Computation Time = %fs\n Error = %f\n', n,dt,tge_iter,ttge,new_error))
saveas(gca,'Transient_GS_ex.jpeg')
end

% Transient_Explicit_SOR on Gauss-Seidel
function T_in = sor_gs_iter_tse(n, new_error, tolerance_error, alpha, dx, dy, xx, yy, omega, T_in, T_old)
dt = 0.01;
k1 = (alpha * dt)/dx^2;
k2 = (alpha * dt)/dy^2;
sge_iter = 1;
nt = 200;

tic;
for tl = 1:nt
    while(new_error > tolerance_error)
        for i = 2:n-1
            for j = 2:n-1
                 T_in(i,j) = ((1-omega)*T_old(i,j)) + (omega*(((1-(2*k1)-(2*k2))*T_old(i,j)) + ((T_in(i-1,j)+T_old(i+1,j))*k1) + ((T_in(i,j-1)+T_old(i,j+1))*k2)));
            end
        end
    new_error = max(max(abs(T_old-T_in)));
    T_old = T_in;
    sge_iter = sge_iter + 1;
    end
end
ttsge = toc;

[c,d] = contourf(xx,yy,T_in);
colormap(jet)
colorbar
xlabel('X-axis')
ylabel('Y-axis')
clabel(c,d);
colorbar
title(sprintf('Transient state (Explicit) using SOR on Gauss-Seidel iterative solver\n Number of nodes = %d  ; Timestep = %.2f\n Number of iteration= %d ; Computation Time = %fs\n Error = %f\n', n,dt,sge_iter,ttsge,new_error))
saveas(gca,'Transient_SOR_GS_ex.jpeg')
end


% Implicit Functions

% Transient_Implicit_Jacobi
function T_in = jacobi_iter_tsi(new_error, tolerance_error, alpha, n, xx, yy, dx, dy, T_in, T_old, T_prev_timestep)
dt = 100;
k1 = (alpha * dt)/dx^2;
k2 = (alpha * dt)/dy^2;
tji_iter = 1;
nt = dt/2;

tic;
for tl = 1:nt
    while(new_error > tolerance_error)
        for i = 2:n-1
            for j = 2:n-1
                T_in(i,j) = (T_prev_timestep(i,j) + (k1*((T_old(i-1,j)) + (T_old(i+1,j)))) + (k2*((T_old(i,j-1)) + (T_old(i,j+1)))))/(1+(2*k1)+(2*k2));
            end
        end
    new_error = max(max(abs(T_old-T_in)));
    T_old = T_in;
    tji_iter = tji_iter + 1;
    end
T_prev_timestep = T_in;
end
ttji = toc;

[c,d] = contourf(xx,yy,T_in);
colormap(jet)
colorbar
xlabel('X-axis')
ylabel('Y-axis')
clabel(c,d);
colorbar
title(sprintf('Transient state (Implicit) using Jacobi iterative solver\n Number of nodes = %d  ; Timestep = %.2f\n Number of iteration= %d ; Computation Time = %fs\n Error = %f\n', n,dt,tji_iter,ttji,new_error))
saveas(gca,'Transient_Jac_im.jpeg')
end

% Transient_Implicit_Gauss Seidel
function T_in = gs_iter_tsi(new_error, tolerance_error, alpha, n, xx, yy, dx, dy, T_in, T_old, T_prev_timestep)
dt = 100;
k1 = (alpha * dt)/dx^2;
k2 = (alpha * dt)/dy^2;
tgi_iter = 1;
nt = dt/2;

tic;
while(new_error > tolerance_error)
for tl = 1:nt
    for i = 2:n-1
        for j = 2:n-1
            T_in(i,j) = (T_prev_timestep(i,j) + (k1*((T_in(i-1,j)) + (T_old(i+1,j)))) + (k2*((T_in(i,j-1)) + (T_old(i,j+1)))))/(1+(2*k1)+(2*k2));
        end
    end
new_error = max(max(abs(T_old-T_in)));
T_old = T_in;
tgi_iter = tgi_iter + 1;
end
T_prev_timestep = T_in;
end
ttgi = toc;

[c,d] = contourf(xx,yy,T_in);
colormap(jet)
colorbar
xlabel('X-axis')
ylabel('Y-axis')
clabel(c,d);
colorbar
title(sprintf('Transient state (Implicit) using Gauss-Seidel iterative solver\n Number of nodes = %d  ; Timestep = %.2f\n Number of iteration= %d ; Computation Time = %fs\n Error = %f\n', n,dt,tgi_iter,ttgi,new_error))
saveas(gca,'Transient_GS_im.jpeg')
end

% Transient_Implicit_SOR on Gauss-Seidel
function T_in = sor_gs_iter_tsi(new_error, tolerance_error, alpha, n, xx, yy, dx, dy, omega, T_in, T_old, T_prev_timestep)
dt = 100;
k1 = (alpha * dt)/dx^2;
k2 = (alpha * dt)/dy^2;
sgi_iter = 1;
nt = dt/2;

tic;
while(new_error > tolerance_error)
for tl = 1:nt
    for i = 2:n-1
        for j = 2:n-1
            T_gsi(i,j) = (T_prev_timestep(i,j) + (k1*((T_in(i-1,j)) + (T_old(i+1,j)))) + (k2*((T_in(i,j-1)) + (T_old(i,j+1)))))/(1+(2*k1)+(2*k2));
            T_in(i,j) = ((1-omega)*T_old(i,j)) + (omega*(T_gsi(i,j)));
        end
    end
new_error = max(max(abs(T_old-T_in)));
T_old = T_in;
sgi_iter = sgi_iter + 1;
end
T_prev_timestep = T_in;
end
ttsgi = toc;

[c,d] = contourf(xx,yy,T_in);
colormap(jet)
colorbar
xlabel('X-axis')
ylabel('Y-axis')
clabel(c,d);
colorbar
title(sprintf('Transient state using (Implicit) SOR on Gauss-Seidel iterative solver\n Number of nodes = %d  ; Timestep = %.2f\n Number of iteration= %d ; Computation Time = %fs\n Error = %.9f\n', n,dt,sgi_iter,ttsgi,new_error))
saveas(gca,'Transient_SOR_GS_im.jpeg')
end

3. Output plots

Steady State

Transient State

Explicit Time Integration Scheme

Implicit Time Integration Scheme

4. Inference

From the above plot the following table can be drawn

(1) It can be inferred from the data obtained for each state type that the Successive Over-Relaxation (SOR) iterative solver projects as the most computaional cost effective iterative solver as compared to the regular Jacobi and Gauss-Seidel iterative solvers.

(2) The steady state 2D heat conduction solved using SOR on Gauss-Seidel iterative solver can be seen to produce the fastest convergence since the steady state doesn't account for time marching, therefore doesn't rely on the timestep whereas the transient state relies heavily on the timesteps. Steady state simulation ignores most of the high order terms dealing with the time, thereby not affecting the steady state result as these terms all go zero at steady state, whereas the transient state simulation tends to capture information at each timestep along with dealing with the cross-terms and high order terms. Hence, all these processes taken by a time marching solution in turn takes comparatively more time to converge than a steady state simulation.

(3) If the rate of convergence is not the concern, we can clearly see that the transient state 2D heat conduction equation solved using implicit time integration scheme with SOR on Gauss-Seidel iterative solver gives the least deviation from the absolute error. The error is so small that neither the table nor the plot could display the exact value upto 8 decimal points. Therefore, it can be deduced that even though implicit schemes taken longer computational time, its results are more accurate when compared to the results obtained using explicit scheme.