Solving the 2D heat conduction Steady and Unsteady state equation using Matlab

Objective:

To solve the 2D heat conduction Steady and Unsteady state equation, by different Iterative methods, and solving the Unsteady state equation both explicitly and Implicitly and simulating them through Matlab and Obtaining their respective graphs.

Procedure:

The Steady state heat conduction equation is:

`∂^2T//∂x^2+∂^2T//∂y^2=0`

Solving the equation:

`T_(i,j)^(n+1)={[T_(i−1,j)+T_(i+1,j)]//Δx2+[T_(i,j−1)+T_(i,j+1)]//Δy2}^n//{2(Δx2+Δy2)//Δx2⋅Δy2}`

Here the equation is solved explicitly,i.e., the RHS values are known and hence it is way too easy to solve the equation.

Program for solving the 2D heat conduction Steady state equation by different Iterative techniques:

clear all
close all
clc

%Inputs
nx=10;
ny=nx;
 x=linspace(0,1,nx);
 y=linspace(0,1,ny);
 dx=x(2)-x(1);
 dy=y(2)-y(1);
 
 %Boundary conditions
 T=ones(nx,ny);
 T(1,:)=900;
 T(end,:)=600;
 T(:,1)=400;
 T(:,end)=800;
 Told=T;
 %Conditions
 tol=1e-4;
 error=9e9;
 k=2*(dx^2+dy^2)/(dx^2*dy^2);
 
 Iterative_solver=3;
 
%Jacobi method
if Iterative_solver==1
    tic
    jacobi_iter=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
            error=max(max(abs(Told-T)));
            Told=T;
            
            jacobi_iter=jacobi_iter+1;
            
            %Plotting
            [x1,y1]=contourf(x,y,T);
            colormap(jet)
          
            clabel(x1,y1);
            xlabel('X axis')
            ylabel('Y axis')
            text_title=sprintf('Temp Jacobi Iteration=%d',jacobi_iter)
            
            title(text_title)
            pause(0.001)
        end
    end 
  Time_jacobi=toc;
end


%Gauss seidel method
if Iterative_solver==2
    gs_iter=1
    tic
    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
            error=max(max(abs(Told-T)))
            Told=T;
            
            gs_iter=gs_iter+1;
            
            %Plotting
            [x1,y1]=contourf(x,y,T);
            colormap(jet)
          
            clabel(x1,y1);
            xlabel('X axis')
            ylabel('Y axis')
            text_title=sprintf('Temp GS Iteration=%d',gs_iter)
            
            title(text_title)
            pause(0.001)
        end
    end 
  Time_GS=toc;
end


%SOR method
if Iterative_solver==3
    tic
    SOR_iter=1
    alpha=1.2
 
    while(error>tol)
        for i=2:nx-1
            for j=2:ny-1
                T(i,j)=alpha*((((T(i-1,j)+Told(i+1,j))/dx^2)+((Told(i,j+1)+T(i,j-1))/dy^2))/k)+Told(i,j)*(1-alpha);
            end
            error=max(max(abs(Told-T)));
            Told=T;
            
            SOR_iter=SOR_iter+1;
            
            %Plotting
            [x1,y1]=contourf(x,y,T);
            colormap(jet)
          
            clabel(x1,y1);
            xlabel('X axis')
            ylabel('Y axis')
            text_title=sprintf('Temp SOR Iteration=%d',SOR_iter);
            
            title(text_title);
            pause(0.00001)
        end
    end 
  Time_SOR=toc;
end

Plots obtained by different Iterative techniques:

As we can observe that the Jacobi method undergoes the higher the no. of iterations than that of the other two iterations. The Successive Over-relaxation method takes lesser time than that of the other two and undergoes lesser Iterations to obtain the lesser error than that of error of tolerance.

The Transient or Unsteady state heat conduction equation is:

`∂T//∂t+α⋅[∂^2T//∂x^2+∂^2T//∂y^2]=0`

Solving the equation:

The Transient state or Unsteady state equation can be solved by both explicitly and implicitly as shown below:

Explicitly:

`T_(i,j)^(n+1)=T_(i,j)+α⋅Δt{[T(i−1,j)−2T_(i,j)+T_(i+1,j)_]//Δx2+[T_(i,j−1)−2T_(i,j)+T_(i,j+1)]//Δy2}^n`

The explicitly solved equation can be solved and simulated without any need of any form Iterative methods, but, the time derivative is involved we are including a time loop in the code to solve it.

Implicitly:

`T_(i,j)^(n+1)=T_(i,j)+α⋅Δt{[T(i−1,j)−2T_(i,j)+T_(i+1,j)_]//Δx2+[T_(i,j−1)−2T_(i,j)+T_(i,j+1)]//Δy2}^n+1`

The Implicitly solved equation undergoes various methods pf Iterative solvers like Jacobi, Gauss seidel and Successive over-relaxation methods within a time loop and the graphs are obtained as follows.

Program for solving the 2D heat conduction Transient or Unsteady state equation by different Iterative techniques:

clear all
close all
clc

%Inputs
l=1;
nx=10;
ny=nx;
x=linspace(0,l,nx);
y=linspace(0,l,ny);
dx=x(2)-x(1);
dy=y(2)-y(1);
dt=1e-4;
t=0.3;
nt=t/dt;
alpha=1.4;
tol=1e-4;
error=1;

T=ones(nx,ny);
T(1,:)=900;
T(ny,:)=600;
T(:,1)=400;
T(:,nx)=800;
Told=T;
Tpre=T;



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


Iterative_method=4
%Explicit method
if Iterative_method==1
 for m=1:nt
   
    for i=2:nx-1
        for j=2:ny-1
            T(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 
    Told=T;
    [c,d]=contourf(x,y,T);
    colormap(jet)
    clabel(c,d)
    xlabel('X axis')
    ylabel('Y axis')
    title('Transient state explicit method')
   
    pause(0.000000000000000000000000000000000001)
 end
end



%Implicit Jacobi method
      jacobi_iter=1
      
    
 if Iterative_method==2
   tic
   %Time loop
 for m=1:nt
        %Convergence loop
          while(error>tol)
           %Nodal loop        
           for i=2:nx-1
              for j=2:ny-1
                  T(i,j)=(Tpre(i,j)+k1*(Told(i-1,j)+Told(i+1,j))+k2*(Told(i,j-1)+Told(i,j+1)))/(1+(2*k1)+(2*k2));
              end
                error=max(max(abs(Told-T)));
          Told=T;
          jacobi_iter=jacobi_iter+1;
      %Plotting
     [c,d]=contourf(x,y,T);
     colormap(jet)
     clabel(c,d)
     xlabel('X axis')
     ylabel('Y axis')
     title_text=sprintf('Jacobi imp Iteration=%d',jacobi_iter);
     title(title_text)
          end
       
     Tpre=T;
    
     pause(0.000000000000000000000000001) 
         end 
     
  end
   GS_imp_time=toc;
 end
 
 
 %Implicit GS method
      GS_iter=1
      
 if Iterative_method==3
   tic
   %Time loop
 for m=1:nt
        %Convergence loop
          while(error>tol)
           %Nodal loop        
           for i=2:nx-1
              for j=2:ny-1
                  T(i,j)=(Tpre(i,j)+k1*(T(i-1,j)+Told(i+1,j))+k2*(T(i,j-1)+Told(i,j+1)))/(1+(2*k1)+(2*k2));
              end
                error=max(max(abs(Told-T)));
          Told=T;
          GS_iter=GS_iter+1;
      %Plotting
     [c,d]=contourf(x,y,T);
     colormap(jet)
     clabel(c,d)
     xlabel('X axis')
     ylabel('Y axis')
     title_text=sprintf('GS imp Iteration=%d',GS_iter);
     title(title_text)
          end
       
     Tpre=T;
    
     pause(0.000000000000000000000000001) 
         end 
     
  end
   GS_imp_time=toc;
 end
 
  
 %Implicit SOR method
      SOR_iter=1
      w=1.4;
 if Iterative_method==4
   tic
   %Time loop
 for m=1:nt
        %Convergence loop
          while(error>tol)
           %Nodal loop        
           for i=2:nx-1
              for j=2:ny-1
                  T(i,j)=((Tpre(i,j)+k1*(T(i-1,j)+Told(i+1,j))+k2*(T(i,j-1)+Told(i,j+1)))/(1+(2*k1)+(2*k2)))*w+(1-w)*Tpre(i,j);
              end
                error=max(max(abs(Told-T)));
          Told=T;
          SOR_iter=SOR_iter+1;
      %Plotting
     [c,d]=contourf(x,y,T);
     colormap(jet)
     clabel(c,d)
     xlabel('X axis')
     ylabel('Y axis')
     title_text=sprintf('SOR imp Iteration=%d',SOR_iter);
     title(title_text)
          end
       
     Tpre=T;
    
     pause(0.000000000000000000000000001) 
         end 
     
  end
   SOR_imp_time=toc;
 end

Plots obtained by different Iterative techniques:

Explicitly:

Implicitly:

As we observe the above graphs, Jacobi method undergoes large no. of Iterations than that of other two methods and Successive over-relaxation method takes lesser time and Iterative steps within Implicit method and hence the Explicitly solved equation gives the fastest simulated graphs among all.

Conclusion:

Thus, both Steady and Unsteady state 2D heat conduction equation is solved by different Iterative techniques and hence simulated using Matlab and the Graphs are obtained.