2D HEAT CONDUCTION - STEADY STATE AND TRANSIENT STATE ANALYSIS

OBJECTIVE:

To simulate the 2D heat conduction across a flat plate by programming in matlab.

THEORY:

The 2D heat conduction simulation has been solved by 3 iterative techniques, mainly:

1.Jacobi method

2.Gauss siedel iterative technique

3.Successive Over Relaxation method

Both the steady state and transient state analysis have been performed in this project.

In the transient state analysis, the heat conduction has been solved by two methods,namely:

1.Implicit method

2.Explicit method

In the implicit method, the variables calculated in the iteration equation at a current time step were with the help of variables whose all values were not known at that time step.

Let us consider the following equation:

In the above equation, suppose only X1 at the k+1 th time step is known, then X2 at k+1 th step remains unknown.Thus this type of equation are known as implicit equation where the variables will get solved in the form of iterations.

In the explicit form of equations, the variables solved at a current time step are from values of variables from previous time steps. Hence the right hand side values are all known.

The transient state equation is:

The steady state heat conduction equation is given as:

The second order derivative terms are elliptic in nature and hence and central order differencing scheme has been followed in those case .

The iterative techniques:

1. JACOBI METHOD:

In Jacobi method,the iterative variables were solved using the values of previous iteration.

2.GAUSS SIEDEL METHOD:

In this method, the iterative variables are calculated with the help of values of both present and previous iteration.The values present for the current iteration are preffered and if they are not available then the previous iteration values are used.

3.SUCCESSIVE OVER RELAXATION METHOD:

This iterative method can be used on both the JACOBI and GAUSS SIEDEL iterative technique.It is like a wrapper which enhances the performance of the iterative solvers like GAUSS SIEDEL or JACOBI.

A scaling factor is used in the SOR method which increases or decreases the aggresiveness of the solver. If the scaling factor is less than 1, then the solver is less aggressive and it is termed as under relaxation. While, if the scaling factor is more than 1, it is considered to be more aggresive and is termed as over relaxation.

THE INITIAL TEMPERATURE PROFILE IS:

THE MATLAB CODE FOR STEADY STATE TEMPERATURE DISTERIBUTION IS:

clear all
close all
clc

%inputs
x=linspace(0,1,10);
y=linspace(0,1,10);
dx=x(2)-x(1);
dy=y(2)-y(1);

%calculating the constants for heat conduction
k=2*(2*dx^2)/(dx^2*dy^2);
k1=1/(k*dx^2);
tol=1e-4;

%initializing the temperature array
xy=zeros(10,10);
xy(:,:)=300;

%boundary conditions
xy(1,:)=600;
xy(2:9,1)=400;
xy(10,:)=900;
xy(1:9,10)=800;

%defining meshgrid
T=xy;

%storing values
Told=T;

%stating the iterative method to follow in the steady state
%iter=1 for JACOBI
%iter=2 for GAUSS SIEDEL
%iter=3 for SOR METHOD

iter=1;

%jacobi method
if (iter==1)
  tic
  iteration=1;
   error=9e9;
  while(error>tol)  
   
      for j=2:9
        for i=2:9
            H=Told(i+1,j)+Told(i-1,j);
            V=Told(i,j+1)+Told(i,j-1); 
            T(i,j)=k1*H+k1*V;
        end
    end
   
    error=max(max(abs(Told-T)));
    Told=T;
    iteration=iteration+1;
   
  end 
   %plotting the temperature distribution
   
    contourf(x,y,T)
    xlabel('x axis')
    ylabel('y axis')
    title_text=sprintf('STEADY STATE-JACOBI Iteration number=%d',iteration);
    title(title_text)
    colorbar
    pause(0.003)
    toc
    et1=toc;
end

%gauss siedel method
if (iter==2)
    iteration=1;
  tic
   error=9e9;
  while(error>tol)  
    
      for j=2:9
        for i=2:9
            H=T(i+1,j)+T(i-1,j);
            V=T(i,j+1)+T(i,j-1); 
            T(i,j)=k1*H+k1*V;
        end
    end
   
    error=max(max(abs(Told-T)));
    Told=T;
    iteration=iteration+1;
    
  end 
    %plotting the temperature distribution
   
    contourf(x,y,T)
    xlabel('x axis')
    ylabel('y axis')
    title_text=sprintf('STEADY STATE-GS Iteration number=%d',iteration);
    title(title_text)
    colorbar
    pause(0.003)
    toc
    et2=toc;
end

 %SOR method
if (iter==3)
 tic
  alpha1=1.1;
  iteration=1;
   error=9e9;
 
   while(error>tol)  
    
       for j=2:9
        for i=2:9
            H=T(i+1,j)+T(i-1,j);
            V=T(i,j+1)+T(i,j-1); 
            T(i,j)=(1-alpha1)*Told(i,j)+alpha1*k1*H+alpha1*k1*V;
        end
       end
       error=max(max(abs(Told-T)));
       Told=T;
       iteration=iteration+1;
       
   end
 
%plotting the temperature distribution
   
       contourf(x,y,T)
       xlabel('x axis')
       ylabel('y axis')
       title_text=sprintf('STEADY STATE-SOR Iteration number=%d',iteration);
       title(title_text)
       colorbar
       pause(0.003)
    toc
    et3=toc;
end

The steady state condition in heat transfer means that any increase of heat rate to the body will not affect the temperature of the body.

THE TEMPERATURE DISTRIBUTION PROFILE AFTER STEADY STATE IS ACHIEVED:

1.JACOBI METHOD:

2.GAUSS SIEDEL METHOD:

3. SUCCESSIVE OVER RELAXATION METHOD:

COMPARING THE FASTEST STEADY STATE SIMULATION:

The simulation time depends on the number of iterations performed by the solver to reach closer to the exact answer.

In the above plot it can be seen that the JACOBI SOLVER has taken the most time to converge whereas the SOR SOLVER has taken the least time to converge.

JACOBI METHOD:208

GS METHOD:112

SOR METHOD:93

THE MATLAB CODE FOR TRANSIENT STATE TEMPERATURE DISTRIBUTION:

IMPLICIT SCHEME:

%transient heat conduction by implicit method

clear all 
close all
clc

%inputs
x=linspace(0,1,10);
y=linspace(0,1,10);
dx=x(2)-x(1);
dy=y(2)-y(1);

%defining constants of heat conduction
alpha=1.2;
alpha1=1.1;
t=1;
C=0.2;

%calculating constants of heat conduction
dt=C*(dx^2)/(2*alpha);
nt=t/dt;
k1=C/2;
k2=C/2;
term1=1/(1+(2*k1)+(2*k2));
term2=k1*term1;
term3=k2*term1;

%initializing the temperature array
xy=zeros(10,10);
xy(:,:)=300;

%defining boundary conditions

xy(1,:)=600;
xy(2:9,1)=400;
xy(10,:)=900;
xy(1:9,10)=800;

%defining the meshgrids
T=xy;

%storing values
Told=T;
Tprev_dt=T;

%defining error and tolerances
tol=1e-4;
error=9e9;

%stating the iterative method to follow in the implicit scheme
%iter=1 for JACOBI
%iter=2 for GAUSS SIEDEL
%iter=3 for SOR METHOD

iter=3;   

%jacobi
  if(iter==1)
      tic
      iteration=1;
     
       %time loop
       
      for k=1:nt
         error=90;
         
          while (error>tol)         
             %space loop
               
              for j=2:9
                  for i=2:9
                    H=Told(i-1,j)+Told(i+1,j);
                    V=Told(i,j+1)+Told(i,j-1);
                    
                    T(i,j)=(term1*Tprev_dt(i,j))+(term2)*H+(term3)*V;
                  end
              end
              
               error=max(max(abs(Told-T)));
               Told=T;
               iteration=iteration+1;               
           end
        Tprev_dt=T;
      end
      contourf(x,y,T)
      xlabel('x axis')
      ylabel('y axis')
      title_text=sprintf('TRANSIENT STATE IMPLICIT-JACOBI METHOD Iteration number=%d',iteration);
      title(title_text)
      colorbar
      
  end
  
  et=toc;
  
  %gauss siedel
   if (iter==2)
       tic
       iteration=1;
       
        %time loop
        
        for k=1:nt
            error=90;
            
             while (error>tol)
                 
               %space loop
               for j=2:9
                  for i=2:9
                    H=T(i-1,j)+T(i+1,j);
                    V=T(i,j+1)+T(i,j-1);
                    T(i,j)=(term1*Tprev_dt(i,j))+(term2)*H+(term3)*V;
                  end
               end
               
               error=abs(max(max(Told-T)));
               Told=T;
               iteration=iteration+1;              
           end
        Tprev_dt=T;
        end
       contourf(x,y,T)
       xlabel('x axis')
       ylabel('y axis')
       title_text=sprintf('TRANSIENT STATE IMPLICIT-GAUSS SIEDEL Iteration number=%d',iteration);
       title(title_text)
       colorbar
       pause(0.003) 
   end
   toc
   et=toc;
   
   
%Succesive over relaxation 
  
 if (iter==3)
     tic
    %time loop
    
    iteration=1;
     for k=1:nt
         error=9e9;
         alpha1=1.1;
         
         while (error>tol)
             
             %space loop
              
             for j=2:9
                  for i=2:9
                    H=T(i-1,j)+T(i+1,j);
                    V=T(i,j+1)+T(i,j-1);
                   T(i,j)=((1-alpha1)*Tprev_dt(i,j))+((alpha1*term1*Tprev_dt(i,j))+(alpha1*((term2)*H+(term3)*V)));
                  end
               end
              
               error=max(max(abs(Told-T)));
               Told=T;
               iteration=iteration+1;
                
           
      
       
           end
        Tprev_dt=T;
     end
     contourf(x,y,T)
     xlabel('x axis')
     ylabel('y axis')
     title_text=sprintf('TRANSIENT STATE IMPLICIT-SOR METHOD Iteration number=%d',iteration);
     title(title_text)
     colorbar
     pause(0.03)
 end 
 et=toc;

The temperature of the different points in the 2D slab keeps changing with time in case of transient state.

THE TEMPERATURE DISTRIBUTION PROFILE IN CASE OF TRANSIENT STATE:

1.JACOBI METHOD:

2.GAUSS SIEDEL METHOD:

3.SOR METHOD:

THE MATLAB CODE FOR TRANSIENT STATE TEMPERATURE DISTRIBUTION:

EXPLICIT SCHEME:

clear all 
close all
clc

%inputs
x=linspace(0,1,10);
y=linspace(0,1,10);
dx=x(2)-x(1);
dy=y(2)-y(1);

%defining constants of heat conduction
alpha=1.2;
t=1;
C=0.2;

%calculating constants of heat conduction
dt=C*(dx^2)/(2*alpha);
nt=t/dt;

%initializing the temperature array
xy=zeros(10,10);
xy(:,:)=300;

%defining boundary conditions

xy(1,:)=600;
xy(2:9,1)=400;
xy(10,:)=900;
xy(1:9,10)=800;

%defining the meshgrids
T=xy;

%storing values
Told=T;

%defining error and tolerances
tol=1e-4;
error=9e9;

%stating the iterative method to follow in the explicit scheme

%jacobi method
iter=1;

%time loop
for k=1:nt
    error=90;
   while (error>tol)
      
       %space loop
       for j=2:9
           for i=2:9
               H=Told(i-1,j)-2*Told(i,j)+Told(i+1,j);
               V=Told(i,j+1)-2*Told(i,j)+Told(i,j-1);
               T(i,j)=Told(i,j)+(C/2)*H+(C/2)*V;
           end
       end
           
       error=max(abs(Told-T));
       Told=T;
       iter=iter+1;
   end
  
end
       contourf(x,y,T)
       xlabel('x axis')
       ylabel('y axis')
       title_text=sprintf('TRANSIENT EXPLICIT-Iteration number=%d',iter);
       title(title_text)
       colorbar
       pause(0.03)

THE TEMPERATURE DISTRIBUTION IN TRANSIENT STATE EXPLICIT SCHEME IS:

COMPARING THE FASTEST TRANSIENT STATE SIMULATION:

It can be seen that JACOBI METHOD has the highest simulation time whereas the EXPLICIT and GS METHOD has the least simulation time.

CHECKING FOR THE STABILITY ANALYSIS:

The stability of a given solution can be analysed by the CFL CONDITION. The COURANT-FRIEDREICHS-LEWY condition is necessary for checking on the convergence of the solutions of partial differential equations.

The COURANTS NUMBER is given as:

IN EXPLICIT METHOD:

When the CFL no is less than 0.5, (here it is taken as 0.2), it is seen that the solution is completely perfect.

But when a CFL number greater than 0.5 is taken,(here it is 0.6), then the solution is seen to be completely blowing up.

IN IMPLICIT METHOD:

In implicit method the CFL number is not seen to have much effect because in both cases of less than or greater than 0.5 CFL number, the solution is found to be stable.

It is therefore can be seen that the implicit scheme is highly stable for any CFL number and steady stae problems can hence for steady state simulation we can easily go for implicit schemes.