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

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

RRaghavendra Kadage8 Feb 2022

** 2D Steady Heat Conduction

*Main code:

clear all
close all
clc
%% Define domain
nx=21;        % Number of points in x direction
ny=21;        % Number of points in y direction
L=1;
H=1;
x=linspace(0,L,nx);
y=linspace(0,H,ny);
dx=L/(nx-1);
dy=H/(ny-1);
M=(dy^2)/(2*(dx^2+dy^2));
N=(dx^2)/(2*(dx^2+dy^2));
omega=1.5;
%% Define Boundary
T=zeros(ny,nx);
T(2:ny-1,1)=400;  %Left
T(2:ny-1,nx)=800; %Right
T(1,2:nx-1)=600;  %Top
T(ny,2:nx-1)=900; % Bottom
Told=T;
%% Inputs
error=1;
tolerance=1e-6;
iteration1=0;
iteration2=0;
iteration3=0;

%% Jacobi
[T_Jacobi]=steady_jacobi(nx,ny,x,y,dx,dy,T,Told,error,tolerance,iteration1,M,N);


%% Gauss Seidal
[T_GS]=steady_gs(nx,ny,x,y,dx,dy,T,Told,error,tolerance,iteration2,M,N);

%% SOR
[T_SOR]=steady_sor(nx,ny,x,y,dx,dy,T,Told,error,tolerance,iteration3,M,N,omega);

1. Jacobi method function:

function[T]=steady_jacobi(nx,ny,x,y,dx,dy,T,Told,error,tolerance,iteration1,M,N)
%% Loops
while error>tolerance      % Loop controller
    for j=2:nx-1
        for i=2:ny-1
            T(i,j)=(M.*(Told(i+1,j)+Told(i-1,j))) + (N.*(Told(i,j+1)+Told(i,j-1)));
        end
    end
        error= max(max(abs(Told-T)));
        Told=T;
        iteration1=iteration1+1;
        


end
        figure(1)
        x_dom=((1:nx)-1).*dx;
        y_dom=(1-((1:ny)-1).*dy);
        [X,Y]=meshgrid(x_dom,y_dom);
        contourf(X,Y,T,25,'ShowText','ON');
        colorbar
        title(sprintf('Steady State Heat Conduction Jacobi Iterations = %d',iteration1))
end

Temperature contour:

2. Gauss Seidal method function:

function[T]=steady_jacobi(nx,ny,x,y,dx,dy,T,Told,error,tolerance,iteration2,M,N)
%% Loops
while error>tolerance
    for j=2:nx-1
        for i=2:ny-1
            T(i,j)=(M.*(T(i+1,j)+T(i-1,j))) + (N.*(T(i,j+1)+T(i,j-1)));
        end
    end
        error= max(max(abs(Told-T)));
        Told=T;
        iteration2=iteration2+1;
        


end
        figure(2)
        x_dom=((1:nx)-1).*dx;
        y_dom=(1-((1:ny)-1).*dy);
        [X,Y]=meshgrid(x_dom,y_dom);
        contourf(X,Y,T,25,'ShowText','ON');
        colorbar
        title(sprintf('Steady State Heat Conduction Gauss Seidal Iterations=%d',iteration2))
end

Temperature contour:

3. SOR function code:

function[T]=steady_sor(nx,ny,x,y,dx,dy,T,Told,error,tolerance,iteration3,M,N,omega)
%% Loops
while error>tolerance
    for j=2:nx-1
        for i=2:ny-1
            TGS(i,j)=(M.*(T(i+1,j)+T(i-1,j))) + (N.*(T(i,j+1)+T(i,j-1)));
            T(i,j)=(1-omega)*Told(i,j)+omega*TGS(i,j);
        end
    end
         error= max(max(abs(Told-T)));
        Told=T;
        iteration3=iteration3+1;
        


end
        figure(3)
        x_dom=((1:nx)-1).*dx;
        y_dom=(1-((1:ny)-1).*dy);
        [X,Y]=meshgrid(x_dom,y_dom);
        contourf(X,Y,T,25,'ShowText','ON');
        colorbar
        title(sprintf('Steady State Heat Conduction SOR Iterations=%d',iteration3))

end

Temperature contour:

* Unsteady Explicit:

Main code:

clc
clear all
close all


%% Create domain
L=1;
B=1;

nx=21;
ny=21;
dx=L/(nx-1);
dy=B/(ny-1);

x=linspace(0,L,nx);
y=linspace(0,B,ny);

%% Creating Mesh
T=650*ones(ny,nx);
T(2:ny-1,1)=400;  % Left
T(2:ny-1,nx)=800; % Right
T(1,2:nx-1)=600;  % Top
T(ny,2:nx-1)=900; % Bottom
T(1,1)=500;       % Top left corner
T(1,nx)=700;      % Top right corner
T(ny,1)=650;      % Bottom left corner
T(ny,nx)=850;     % Bottom right corner
Tprev=T;

%% Time step
dt=0.001;
%% Thermal diffusivity
alpha=0.5;
%% Control parameters
convergence=1;
tolerance=1e-4;
iteration=0;

%% Special constants
k1=alpha*dt/(dx^2);
k2=alpha*dt/(dy^2);
CFL=k1+k2;

%% Loop
[T]=explicit(nx,ny,x,y,dx,dy,T,Tprev,convergence,tolerance,iteration,k1,k2)

function:

function[T]=explicit(nx,ny,x,y,dx,dy,T,Tprev,convergence,tolerance,iteration,k1,k2)
%% Loops
while convergence>tolerance
    for i=2:ny-1
        for j=2:nx-1
            T(i,j)=(1-2*k1-2*k2)*Tprev(i,j)+k1*(Tprev(i,j+1)+Tprev(i,j-1))+k2*(Tprev(i+1,j)+Tprev(i-1,j));
        end
    end
    convergence=max(max(abs(T-Tprev)));
    Tprev=T;
    iteration=iteration+1;
    
    %% Plotting
        figure(1)
        x_dom=((1:nx)-1).*dx;
        y_dom=(1-((1:ny)-1).*dy);
        [X,Y]=meshgrid(x_dom,y_dom);
        contourf(X,Y,T,25,'ShowText','ON');
        title(sprintf('Explicit Time Iterations=%d',iteration));
end
end

Temperature contour:

* Unsteady implicit:

main code:

clc
clear all
close all
%% Create domain
L=1;
B=1;

nx=11;
ny=11;
dx=L/(nx-1);
dy=B/(ny-1);

x=linspace(0,L,nx);
y=linspace(0,B,ny);

%% Creating Mesh
T=650*ones(ny,nx);
T(2:ny-1,1)=400;
T(2:ny-1,nx)=800;
T(1,2:nx-1)=600;
T(ny,2:nx-1)=900;
T(1,1)=500;
T(1,nx)=700;
T(ny,1)=650;
T(ny,nx)=850;
Told=T;
Tprev=T;
T_initial=T;

%% Time step
dt=0.001;
%% Thermal diffusivity
alpha=0.5;
%% Control parameters
error=1;
convergence=1;
tolerance=1e-4;
space_iteration=0;
time_iteration=0;

%% Special constants
k1=alpha*dt/(dx^2);
k2=alpha*dt/(dy^2);
omega=1.05;

%% Jacobi
[T_jacobi]=jacobi_implicit(nx,ny,x,y,dx,dy,T,Told,Tprev,error,convergence,tolerance,space_iteration,time_iteration,k1,k2);

%% Gauss Seidal
[T_GS]=GS_implicit(nx,ny,x,y,dx,dy,T,Told,Tprev,error,convergence,tolerance,space_iteration,time_iteration,k1,k2);

%% SOR
[T_SOR]=SOR_implicit(nx,ny,x,y,dx,dy,T,Told,Tprev,error,convergence,tolerance,space_iteration,time_iteration,k1,k2,omega);

1. Jacobi method function code:

function[T_Jacobi]=jacobi_implicit(nx,ny,x,y,dx,dy,T,Told,Tprev,error,convergence,tolerance,space_iteration,time_iteration,k1,k2)

%% Loops
% time loops
while convergence>tolerance
    % space loops
    while error>tolerance
        for i=2:ny-1
            for j=2:nx-1
                T(i,j)=(1+2*k1+2*k2)^-1*(Tprev(i,j)+(k1*(Told(i,j+1)+Told(i,j-1)))+(k2*(Told(i+1,j)+Told(i-1,j))));
            end
        end
        error=max(max(abs(T-Told)));
        Told=T;
        space_iteration=space_iteration+1;
    end    
    convergence=max(max(abs(T-Tprev)));
    Tprev=T;
    time_iteration=time_iteration+1;
    error=1;
    iteration=0;
        
end
        figure(1)
        x_dom=((1:nx)-1).*dx;
        y_dom=(1-((1:ny)-1).*dy);
        [X,Y]=meshgrid(x_dom,y_dom);
        contourf(X,Y,T,25,'ShowText','ON');
        colorbar;
        title(sprintf('Jacobi Time Iterations=%d \n Jacobi Space Iterations =%d',time_iteration,space_iteration))
        T_Jacobi=T;
end

Temperature contour:

2. Gauss Sedial method function code:

function[T_GS]=GS_implicit(nx,ny,x,y,dx,dy,T,Told,Tprev,error,convergence,tolerance,space_iteration,time_iteration,k1,k2,omega)

%% Loops
% time loops
while convergence>tolerance
    % space loops
    while error>tolerance
        for i=2:ny-1
            for j=2:nx-1
                T(i,j)=(1+2*k1+2*k2)^-1*(Tprev(i,j)+(k1*(T(i,j+1)+T(i,j-1)))+(k2*(T(i+1,j)+T(i-1,j))));
               
            end
        end
        error=max(max(abs(T-Told)));
        Told=T;
        space_iteration=space_iteration+1;
    end    
    convergence=max(max(abs(T-Tprev)));
    Tprev=T;
    time_iteration=time_iteration+1;
    error=1;
    iteration=0;
end
        figure(2)
        x_dom=((1:nx)-1).*dx;
        y_dom=(1-((1:ny)-1).*dy);
        [X,Y]=meshgrid(x_dom,y_dom);
        contourf(X,Y,T,25,'ShowText','ON');
        colorbar;
        title(sprintf('GS Time Iterations=%d \n Space Iterations=%d',time_iteration,space_iteration))
        T_GS=T;
        
end

Temperature control:

3. SOR method function code:

function[T_SOR]=SOR_implicit(nx,ny,x,y,dx,dy,T,Told,Tprev,error,convergence,tolerance,space_iteration,time_iteration,k1,k2,omega)

%% Loops
% time loops
while convergence>tolerance
    % space loops
    while error>tolerance
        for i=2:ny-1
            for j=2:nx-1
                T_GS(i,j)=(1+2*k1+2*k2)^-1*(Tprev(i,j)+(k1*(T(i,j+1)+T(i,j-1)))+(k2*(T(i+1,j)+T(i-1,j))));
                T(i,j)= (1-omega)*Told(i,j) + omega*T_GS(i,j);
            end
        end
        error=max(max(abs(T-Told)));
        Told=T;
        space_iteration=space_iteration+1;
    end    
    convergence=max(max(abs(T-Tprev)));
    Tprev=T;
    time_iteration=time_iteration+1;
    error=1;
    iteration=0;
        
end
        figure(3)
        x_dom=((1:nx)-1).*dx;
        y_dom=(1-((1:ny)-1).*dy);
        [X,Y]=meshgrid(x_dom,y_dom);
        contourf(X,Y,T,25,'ShowText','ON');
        colorbar;
        title(sprintf('SOR Time Iterations=%d \n Space Iterations=%d',time_iteration,space_iteration))
        T_SOR=T;
end

Temperature contour

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.