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

2D Conduction Solution-Steady and Transient Analysis

DDinesh Hase14 Nov 2019

-Steady and Unsteady 2D Heat conduction solution for unit square domain (x=y)

  • 2-D Steady Heat Conduction equation

`(del^2T)/(del^2x) +(del^2T)/(del^2y) =0`

  • 2-D UnSteady Heat Conduction equation

`(delT)/(delt)+alpha((del^2T)/(del^2x) +(del^2T)/(del^2y) )=0`

where alpha = thermal diffusivity of material

t = time

T =Temperature

1) 2-D Steady Heat Conduction equation

`(del^2T)/(del^2x) +(del^2T)/(del^2y) =0`

solving above equation using central difference scheme for both x and y direction.

Following 3 Iterative solver methods are used.

  • Jacobi Method
  • Gauss-Siedel Method
  • Successive over relaxation Method (SOR)

Code for steady State Condition

%Steady state analysis of 2-D conduction equation using 3 iterative technique     
clear all
close all
clc
%Defining system
n_x=100;
n_y=n_x;%square domain
l=1;%length of domain
x=linspace(0,l,n_x);
%y=linspace(0,l,n_y);
y=x;
[xx,yy]=meshgrid(x,y);
dx=x(2)-x(1);
dy=dx;
%defining B.C.
T=298*ones(n_x,n_y);
T(:,1)=400;%in kelvin
T(1,:)=600;
T(:,n_y)=800;
T(n_x,:)=900;
T_old=T;
error=9e9;%for first iteration
error_limit=1e-4;
iteration_method=3;%Select method of iteration

if (iteration_method==1)
%Jacobi Method
jacobi_iteration =1;
tic
error=9e9
while error > error_limit
    
for i=2:n_x-1
    
    for j=2:n_y-1
        
        Tl=T_old(i,j-1);
        Tr=T_old(i,j+1);
        Tt=T_old(i+1,j);
        Tb=T_old(i-1,j);
        
        T(i,j)=((Tr+Tl)+(Tb+Tt))*0.25;
    end   
end
jacobi_iteration=1+jacobi_iteration;
error=max(max(abs(T_old-T)));
T_old=T;
end
Time_of_iteration=toc;
[c,h]=contourf(xx,yy,T,\'--\')
set(h,\'ShowText\',\'on\')
colorbar
colormap(\'jet\')
xlabel(\'x coordinate\')
ylabel(\'y coordinate\')
title(sprintf(\'2-D Steady Conduction jocobi Method \\n no. of iteration=%d \\n Time=%0.2fsec\',jacobi_iteration,Time_of_iteration))
end

%%%%%%%%%%%
if (iteration_method==2)
%Gauss Siedel Method
gs_iteration =1;
tic
while error > error_limit
    
for i=2:n_x-1
    
    for j=2:n_y-1
        
        Tl=T(i,j-1);
        Tr=T(i,j+1);
        Tt=T(i+1,j);
        Tb=T(i-1,j);
        
        T(i,j)=((Tr+Tl)+(Tb+Tt))*0.25;
    end   
end
gs_iteration=1+gs_iteration;
error=max(max(abs(T_old-T)));
T_old=T;
end
Time_of_iteration=toc;
[c,h]=contourf(xx,yy,T,\'--\')
set(h,\'ShowText\',\'on\')
colorbar
colormap(\'jet\')
xlabel(\'x coordinate\')
ylabel(\'y coordinate\')
title(sprintf(\'2-D Steady Conduction Gauss Siedel Method \\n no. of iteration=%d\\n Time=%0.2fsec\',gs_iteration,Time_of_iteration))
end

%%%%%%%%%%%%%
if (iteration_method==3)
%SOR Method
sor_iteration =1;
omega=0.9;%under relaxation
tic
while error > error_limit
    
for i=2:n_x-1
    
    for j=2:n_y-1
        
        Tl=T(i,j-1);
        Tr=T(i,j+1);
        Tt=T(i+1,j);
        Tb=T(i-1,j);
        
        T(i,j)=(T_old(i,j)*(1-omega))+omega*(((Tr+Tl)+(Tb+Tt))*0.25);
        %in above expression GS method is used to find slope hencevif we
        %use alpha=1 then no of iteration are same as gs iteration.
    end   
end
sor_iteration=1+sor_iteration;
error=max(max(abs(T_old-T)));
T_old=T;
end
Time_of_iteration=toc;
[c,h]=contourf(xx,yy,T,\'--\')
set(h,\'ShowText\',\'on\')
colorbar
colormap(\'jet\')
xlabel(\'x coordinate\')
ylabel(\'y coordinate\')
title(sprintf(\'2-D Steady Conduction SOR Method \\n no. of iteration=%d \\n Time=%0.2fsec\',sor_iteration,Time_of_iteration))
end

Steady State plots

1)Jacobi Method

2)Gauss Siedel Method

3)SOR Method

omega=1.9

omega=0.9

As per the plots above the no. of iteration are least for SOR method and decreases with increase in omega.

For jacobi, for each iteration we are using the predefined value of temperature for each node while

for Gauss-Siedel method we are using the calculated value of temperature for next iteration.

-Plot without convergence loop(while loop)

Steady state means already temperature is in equilibrium i.e. heat transfer in = heat transfer out.hence as we change the error limit we will get different plots.

2) 2-D Transient (UnSteady) Heat Conduction equation =

`(delT)/(delt)+alpha((del^2T)/(del^2x) +(del^2T)/(del^2y) )=0`

Transient= Temperature changing w. r. t. time

Transient analysis can be done by using

  • Explicit Method
  • Implicit Method

1) Explicit Method

In Explicit method we have equation for each iteration with only one unknown. For that we don\'t required convergence loop.

In Explicit Transient analysis we are calculating temperature at next time step from the previous time step.

`Tp^(n+1)=Tp^n+(alpha*dt/dx^2 *(Tr+Tl+Tt+Tb+Tp)^n)`

%Transient analysis of 2-D conduction equation using 
iterative technique-Explicit
clear all
close all
clc

%Defining system
n_x=100;
n_y=n_x;
l=1;
x=linspace(0,l,n_x);
%y=linspace(0,l,n_y);
y=x;
[xx,yy]=meshgrid(x,y);
dx=x(2)-x(1);
dy=dx;
%defining B.C.
T=298*ones(n_x,n_y); 
T(:,1)=400;
T(1,:)=600;
T(:,n_y)=800;
T(n_x,:)=900;
T_old=T;
T_const=T_old;
error=9e9;
error_limit=1e-4;
alpha=0.000127;%thermal diffusivity of gold
t=1000;
dt=0.1;
n_t=t/dt;
k1=(alpha*dt)/dx^2;
k2=1-(4*k1);
CFL_Number=k1;
iteration_method=3;


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if (iteration_method==1)
%Jacobi Method
jacobi_iteration =1;
tic
for m=1:n_t %time loop
  
for i=2:n_x-1
       for j=2:n_y-1
        
        Tl=T_old(i,j-1);
        Tr=T_old(i,j+1);
        Tt=T_old(i+1,j);
        Tb=T_old(i-1,j);
        
        T(i,j)=(T_old(i,j)*k2)+(k1*(Tr+Tl+Tb+Tt));
       end   
end

jacobi_iteration=1+jacobi_iteration;
T_old=T;
end
Time_of_iteration=toc;
[c,h]=contourf(xx,yy,T,\'--\')
set(h,\'ShowText\',\'on\')
colorbar
colormap(\'jet\')
xlabel(\'x coordinate\')
ylabel(\'y coordinate\')
title(sprintf(\'2-D Transient Conduction Explicit Approach \\n no. of iteration=%d \\n cfl number=%0.2f Time=%0.2fsec\',jacobi_iteration,CFL_Number,Time_of_iteration))
end


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if (iteration_method==2)
%gauss siedel method
gs_iteration =1;
tic
for m=1:n_t %time loop
   
    for i=2:n_x-1
    
    for j=2:n_y-1
        
        Tl=T(i,j-1);
        Tr=T(i,j+1);
        Tt=T(i+1,j);
        Tb=T(i-1,j);
        
        T(i,j)=(T_old(i,j)*k2)+(k1*(Tr+Tl+Tb+Tt));
    end   
end

gs_iteration=1+gs_iteration;
T_old=T;
end
Time_of_iteration=toc;
[c,h]=contourf(xx,yy,T,\'--\')
set(h,\'ShowText\',\'on\')
colorbar
colormap(\'jet\')
xlabel(\'x coordinate\')
ylabel(\'y coordinate\')
title(sprintf(\'2-D Transient Conduction Explicit Approach \\n no. of iteration=%d \\n cfl number=%0.4f Time=%0.2f\',gs_iteration,CFL_Number,Time_of_iteration))
end

if (iteration_method==3)
%SOR method
sor_iteration =1;
alpha=1.1
for k=1:n_t %time loop
  for i=2:n_x-1
    
    for j=2:n_y-1
        
        Tl=T(i,j-1);
        Tr=T(i,j+1);
        Tt=T(i+1,j);
        Tb=T(i-1,j);
        
        T(i,j)=(T_old(i,j)*(1-alpha))+alpha*(T_old(i,j)*k2)+(k1*(Tr+Tl+Tb+Tt));
    end   
end
sor_iteration=1+sor_iteration;
T_old=T;
end

Time_of_iteration=toc;
[c,h]=contourf(xx,yy,T,\'--\');
set(h,\'ShowText\',\'on\')
colorbar
colormap(\'jet\')
xlabel(\'x coordinate\')
ylabel(\'y coordinate\')
title(sprintf(\'2-D Transient Conduction Explicit Approach \\n no. of iteration=%d \\n cfl number=%0.4f Time=%0.4f\',sor_iteration,CFL_Number,Time_of_iteration))
end

Result Transient Explicit-

Result plot for explict approch is same for iterative method

2) Implicit Approach

In implicit Transient analysis in each iteration we are not getting the exact value of temperature rather we are getting equation with more than one unknows.

so we get the equations for each node. So these no.of equaation are solved using different iterative techniques.

`Tp^(n+1)=Tp^n+(alpha*dt/dx^2 *(Tr+Tl+Tt+Tb+Tp)^(n+1))/(1+4*alpha*dt/dx^2)`

%Transient analysis of 2D conduction equation using 3 iterative technique-Implicit
clear all
close all
clc
%Defining system
n_x=100;
n_y=n_x;
l=1;
x=linspace(0,l,n_x);
%y=linspace(0,l,n_y);
y=x;
[xx,yy]=meshgrid(x,y);
dx=x(2)-x(1);
dy=dx;
%defining B.C.
T=298*ones(n_x,n_y);
T(:,1)=400;
T(1,:)=600;
T(:,n_y)=800;
T(n_x,:)=900;
T_old=T;
T_prev=T_old;
T_const=T_old;
error=9e9;
error_limit=1e-4;
alpha=0.000167;
t=1000;
dt=0.1;
n_t=t/dt;
k1=(alpha*dt)/dx^2;
k=1/(1+4*k1);
iteration_method=3;

if (iteration_method==1)
    
%Jacobi Method
jacobi_iteration =1;
tic
for m=1:n_t %time loop
    error=9e9;
while error > error_limit

for i=2:n_x-1
    
    for j=2:n_y-1
        
        Tl=T_old(i,j-1);
        Tr=T_old(i,j+1);
        Tt=T_old(i+1,j);
        Tb=T_old(i-1,j);
        
        T(i,j)=(T_prev(i,j)+k1*(Tr+Tl+Tb+Tt))*k;
    end   
end
jacobi_iteration=1+jacobi_iteration;
error=max(max(abs(T_old-T)));
T_old=T;
end
T_prev =T;
end
Time_of_iteration=toc;
[c,h]=contourf(xx,yy,T,\'--\')
set(h,\'ShowText\',\'on\')
colorbar
colormap(\'jet\')
xlabel(\'x coordinate\')
ylabel(\'y coordinate\')
title(sprintf(\'2-D Transient Conduction Implicit Approach \\n no. of jacobi iteration=%d \\n Time=%0.2fsec\',jacobi_iteration,Time_of_iteration))
end

%%%%%%%%%%
if (iteration_method==2)
%gauss siedel method
gs_iteration=1;
tic
for m=1:n_t %time loop
   error=9e9;
while error > error_limit

for i=2:n_x-1
    
    for j=2:n_y-1
        
        Tl=T(i,j-1);
        Tr=T(i,j+1);
        Tt=T(i+1,j);
        Tb=T(i-1,j);
        
        T(i,j)=(T_prev(i,j)+k1*(Tr+Tl+Tb+Tt))*k;
    end   
end
gs_iteration=1+gs_iteration;
error=max(max(abs(T_old-T)));
T_old=T;
end
T_prev =T;
end
Time_of_iteration=toc
[c,h]=contourf(xx,yy,T,\'--\')
set(h,\'ShowText\',\'on\')
colorbar
colormap(\'jet\')
xlabel(\'x coordinate\')
ylabel(\'y coordinate\')
title(sprintf(\'2-D Transient Conduction Implicit Approach \\n no. of Gauss-Siedel iteration=%d \\n Time=%0.2fsec\',gs_iteration,Time_of_iteration))
end


if (iteration_method==3)
%SOR method 
omega=linspace(0.8,1.2,5);%consider range of omega(0.8,0.9,1,1.1,1.2)

for z=1:length(omega)
    tic
iteration=1; 
for m=1:n_t %time loop
  error=9e9;
 
while error > error_limit
for i=2:n_x-1
    
    for j=2:n_y-1
        
        Tl=T(i,j-1);
        Tr=T(i,j+1);
        Tt=T(i+1,j);
        Tb=T(i-1,j);
        T(i,j)=(T_prev(i,j)*(1-omega(z)))+(omega(z))*(k)*(T_prev(i,j)+k1*(Tr+Tl+Tb+Tt));
    end
end

iteration=1+iteration;
error=max(max(abs(T_old-T)));
T_old=T;
end
T_prev =T;
end
figure(z);
sor_iteration(z)=iteration;
Time_of_iteration(z)=toc;
[c,h(z)]=contourf(xx,yy,T,\'--\');
set(h(z),\'ShowText\',\'on\')
colorbar
colormap(\'jet\')
xlabel(\'x coordinate\')
ylabel(\'y coordinate\')
title(sprintf(\'2-D Transient Conduction Implicit Approach \\n no. of sor iteration=%d \\n time=%0.2f omega=%0.2f\',sor_iteration(z),Time_of_iteration(z),omega(z)))
T_prev=T_const;
end
bar(Time_of_iteration,0.2,\'green\')
set(gca,\'XTicklabel\',{0.8,0.9,1,1.1,1.2})
end

Result for Transient Implicit

1) Jacobi Method

2)Gauss Siedel

3) SOR

In SOR method

omega=0.8

omega=0.9

omega=1.0

omega=1.1

omega=1.2

From above plots the no.of iteration are more for Implicit method than explicit .

In above case, for given input n_x=100, t=1000,dt=0.1,alpha=0.000127, the no,of iteration increases with increse value of omega.

No.of Iteration for different value of omega

Time of Iteration for different value of omega

1)n_x= 100

2)1)n_x= 10

from above bar charts it is clear that optimum value of omega is also function of n_x.

Similarly the optimum value of omega is also depend on the type of problem i.e. 2D or 1D or 3D and domain.For larger value of omega we may get infinite value for temperature i.e. solution is unstable.

when in SOR, omega=1, we are getting same plot with same no. of itearation as of gauss siedel method.

If omega<1 then it is called under relaxation

if omega>1 then it is called as over relaxation.

Also as we increase value of alpha i.e diffusivity then the temperation diffusion is becoming fast.

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.