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

The 2Dheat conduction equation is of 2 types:

(i) steady-state

(ii) unsteady-state

The steady-state equation can be solved with the implicit method while coming to the unsteady state equation there is a time term, so it can be solved either explicit or implicit method.

For steady-state eq

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

where T is temperature and x,y are special coordinates.

`T_p= 1/k*((T_l+T_r)/(Deltax^2)) + 1/k*((T_t+T_b)/(Deltay^2))`

If the above equation in nodes(i,j) form, then the equation should be

`T(i,j) = 1/k*((T_(i-1,j)+T_(i+1,j))/(Deltax^2)) + 1/k*((T_(i,j-1)+T_(i,j+1))/(Deltay^2))`

codes for steady-state equation: 

function code:

 for jacobi

function out = jacobi(t,told,i,j)
   out = 0.25*((told(i,j-1)+told(i,j+1)+told(i-1,j)+told(i+1,j)))
 end

for GS

function out = GS(t,told,i,j)
   out = 0.25*((t(i,j-1)+told(i,j+1)+t(i-1,j)+told(i+1,j)))
end

for Sor

function out = sor(t,told,i,j,f)
  out = (1-f)*told(i,j) + f*(0.25*(t(i,j-1)+told(i,j+1)+t(i-1,j)+told(i+1,j)))
end

 main code:

clear all
close all
clc

% inputs
% number of grid points along x and y axis
nx= 10; 
ny=10;
x=linspace(0,1,nx);
dx = x(2)-x(1);
y= linspace(0,1,ny);
z=fliplr(y);
f=1.1;
iter = 1;

t= 300*ones([nx,ny]);

% applying BC

t(1,:)=600;
t(end,:)=900;
t(:,1)=400;
t(:,end)=800;

% creating a copy
told=t;
r1=t;
 
error= 19;
tol = 1e-4;
for iterative_solver = 1:3;
   while(error>tol);
  for i=2:nx-1;
  for j=2:ny-1;
    if(iterative_solver==1);
      t(i,j)= jacobi(t,told,i,j);    
      elseif(iterative_solver==2)
      t(i,j)=GS(t,told,i,j);
    elseif(iterative_solver==3)
      t(i,j)= sor(t,told,i,j,f);
   end
   end
   end
  error = max(max(t-told));
     told = t;
    [u1,v1]=meshgrid(x,z);
    % plot
    
  if(iterative_solver==1);
  figure(1);
  contourf(u1,v1,t);
  [C,h]= contourf(u1,v1,t); 
  colorbar;
  colormap(jet);
  clabel(C,h);
  xlabel('x-axis');
  ylabel('y-axis');
  title_text=sprintf('jacobi iterantion number=%d',iter);
  title(title_text);
  pause(0.003);
  elseif(iterative_solver==2);
  
  figure(2);
  contourf(u1,v1,t);
  [C,h]= contourf(u1,v1,t); 
  colorbar;
  colormap(jet);
  clabel(C,h);
  xlabel('x-axis');
  ylabel('y-axis');
  title_text=sprintf('GS iterantion number=%d',iter);
  title(title_text);
  pause(0.003);
  else  
  
  figure(3);
  [C,h]= contourf(u1,v1,t); 
  colorbar;
  colormap(jet);
  clabel(C,h);
  xlabel('x-axis');
  ylabel('y-axis');
  title_text=sprintf('sor iterantion number=%d',iter);
  title(title_text);
  pause(0.003);  
  end
  iter=iter+1;
  end
error=19;
iter=1;
t=r1;
told=t;
end

Result:

for jacobi method:

for GS method:

for SOR method:

From the above results, we can conclude that sor method is the fastest method as it solves the equation in 92 iterations.

For unsteady or transient state eq: 

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

where T is time and x,y are special coordinates.

As mentioned before, the unsteady state equation solved in 2 methods. They are Explicit and  Implicit methods.

For Explicit method:  

The equation in terms of nodes(i,j)

`T_p^(n+1)= T_p^n+k1*(T_l-2*T_p+T_r)^n + k2*(T_t-2*T_p+T_b)^n`

where `K_1 = alpha*(trianglet)/(trianglex^2)`,  `K_2 = alpha*(trianglet)/(triangley^2)`

If `K_1 + K_2` is greater than 0.5 then it is a unsteable.

For Implicit method:  

The equation in terms of nodes(i,j)

`T_p^(n+1)= Term1*(T_p^n+k_1*(T_l+T_r)^(n+1) + k_2*(T_t+T_b)^(n+1))`

 where term1 = 1/(1+2*k_1+2*k_2), term2 = `K_1`*term1, term3 = `K_2`*term1

codes for Transient method:

For Explicit method:

clear all
close all
clc

% inputs
% number of grid points along x and y axis
nx= 10; 
ny=10;
x=linspace(0,1,nx);
dx = x(2)-x(1);
y= linspace(0,1,ny);
dy = y(2)-y(1);

c= 1;
% number of steps 
nt=200;
% time steps
dt = 0.001;
% CFL number
k= (c*dt)/dx^2;

t= ones([nx,ny])*673;
% applying BC

t(1,:)=600;
t(end,:)=900;
t(:,1)=400;
t(:,end)=800;

told = t; 

tic
iter=1
for p=1:nt
  for i= 2:nx-1
    for j=2:ny-1
      t(i,j)= k*(told(i-1,j)+told(i+1,j)+told(i,j+1)+told(i,j-1))+ (1-4*k)*told(i,j);
    endfor
  endfor
  error=max(abs(told-t));
  told=t;
  iter=iter+1;
endfor
time = toc;
figure(1)
[c,h]=contourf(x,fliplr(y),t);
colorbar;
colormap(jet);
clabel(c,h);
title(sprintf('2D unsteady heat conduction (Explicit method) n iteration number=%d; t computational time :%0.4f',iter,time));
xlabel('xx');
ylabel('yy');

 Plot:

For Implicit method: 

clear all
close all
clc

% inputs
% number of grid points along x and y axis
nx= 10; 
ny=10;
x=linspace(0,1,nx);
dx = x(2)-x(1);
y= linspace(0,1,ny);
dy = y(2)-y(1);
f= 1.2;

t= ones(nx,ny)*298;
% applying BC

t(1,:)=600;
t(end,:)=900;
t(:,1)=400;
t(:,end)=800;

alpha= 1.4;
dt=0.001;
 
% defining k1 and k2
k1=alpha*(dt/dx^2);
k2=alpha*(dt/dy^2);

told=t;
t_initial=t;

for iterartive_solver = 1:3
 nt=1:1400;
 error=9e9;
 tol=1;
 tic;
 
% iterartive_solver==1:3;
  iter=1;
 while error>tol
   for i=2:nx-1
     for j=2:ny-1
       term1= 1/(1+2*k1+2*k2);
       term2= k1*term1;
       term3= k2*term1; 
       if iterartive_solver==1;
       t(i,j)= ((told(i,j)*term1)+(term2*(t(i-1,j)+t(i+1,j))+(term3*(t(i,j-1)+t(i,j+1)))));
     elseif iterartive_solver==2;
       t(i,j)= ((told(i,j)*term1)+(term2*(t(i-1,j)+told(i+1,j))+(term3*(t(i,j-1)+told(i,j+1)))));
      elseif iterartive_solver==3;
       t(i,j)= (1-f)*told(i,j)+ f*(((told(i,j)*term1)+(term2*(t(i-1,j)+told(i+1,j))+(term3*(t(i,j-1)+told(i,j+1))))));
     end
   end
   error_n=max(abs(told-t));
   error=max(error_n);
   told=t;
   iter=iter+1;
   
 end
end
time=toc;

%plot

if iterartive_solver==1;
figure(1)
contourf(x,fliplr(y),t);
[C,h]=contourf(x,fliplr(y),t);
colorbar;
colormap(jet);
clabel(C,h);
text = sprintf('JACOBI Method n iteration number = %d; t computational time: %0.4f',iter ,time);
title(text);
told=t_initial;
t=t_initial;
elseif iterartive_solver==2;

figure(2)
[c,h]=contourf(x,fliplr(y),t);
colorbar;
colormap(jet);
clabel(c,h);
title(sprintf('GS Method n iteration number = %d; t computational time: %0.4f',iter,time));
told=t_initial;
t=t_initial;
else

figure(3)
[c,h]=contourf(x,fliplr(y),t);
colorbar;
colormap(jet);
clabel(c,h);
title(sprintf('SOR Method n iteration number = %d; t computational time: %0.4f',iter,time));
told=t_initial;
t=t_initial;
end 
end

Plots:

for Jacobi:

for GS:

for SOR: 

Explanation: 

While comparing both Explicit and Implicit methods results, the Explicit method is much quicker than the implicit method.

The Explicit method takes 201 iterations whereas the Implicit method takes 577, 577, 521 iterations.

Coming to the implicit method SOR is the fastest method while comparing it with Jacobi and GS. The SOR method took 521 iterations.

But the Implicit method is more accurate as it uses various iterative solvers.