Challenge

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

MMandalemula Sai Baba14 Aug 2022

Aim:- Solving 2D Transient heta conduction equation.

# For Un-steady 2D heat conduction equation is reprented as ..... `del`T/`del`t = `alpha`(`del`^2T/`del`x^2+`del`^2T/`del`y^2)...........

# for steady state equation, the heat conduction wquation is reprensented by

# `del`^2/`del`x^2 + `del`^2T/`del`y^2)..............

# These equations are partially differential equation and hence the numerical discretizations equations can be used as numerical analysis.

# Boundary conditions

# Top boundary = 600k

# right boundary = 800k

# left boundary = 400k

# bottom boundary = 900k

# For transient analysis , either the equations can be solved in its implicit and explicit form.

# `del`T/`del`t is discretized using forward time maching scheme as Ti,j^n+1-Ti,j^n/`delta`t...............................

# (i,j) denotes the nodes and the super script of the time step.

# since the right hand side of the equation is a second order partially derivative, the central differencing scheme is used for discretizing interms as follows.

# `del`^2.T/`del`x^2 = Ti-1,j-2Ti,j+Ti+1,j/`delta`x^2................

# `delta`^2T/`delta`y^2 = Ti,j-1-2Ti,j+Ti,j+1/`delta`y^2....................................

# Depending on which tie step these second order derivative are computed , they can be termed as implicit and explict form...

# EXPLICIT FORM:-

# When computed at the previous time in the form equation is the explict form

# T_p^n+1-T_p^n/`delta`t = `alpha`[T_L^n-2T_p^n+T_R^n/dx^2 + T_T^n-2T_P^n+T_B^n/dy^2]..........................

# rearriging:- T_p^n+1 = T_p^n+`alpha`(dt) .[[T_L^n-2T_p^n+T_R^n/dx^2 + T_T^n-2T_P^n+T_B^n/dy^2]......................

let us consider

# alpha*dt/dy^2 = k1

# alpha *dt/dy^2 = K2...............

# substuting the values of given equation

# T_p^n+1 = T_p^n+k1*(T_L^n-2T_p^n+T_R^n)+k2*(T_T^n-2T_p^n+T_b^n)................................

#T_p^n+1 = T_P^n(1-2k1+k2)+k1(T_L^n+T_R^n)+k2(T_T^n+T_B^n)...........................................................

IMPLICIT FORM:-

T_P^n+1-T_P^n/dt = alpha*[T_L^n+1-2.T_P^n+1+T_R^n+1/dx^2+T_T^n+1-2.T_p^n+1+T_B^n+1/dy^2].....................................

# T_P^n+1 = T_P^n+alpha(dt) [T_L^n+1-2.T_P^n+1+T_R^n+1/dx^2+T_T^n+1-2.T_p^n+1+T_B^n+1/dy^2]........................

Rearranging:- T_p^n+1 = T_p^n+k1*(T_L^n+1+T_R^n+1)+k2(T_T^n+1+T_B^n+1)/(1+2.k1*k2)............................................................

steady state:-

jacobi method octave code:-

clear all
close all
domain_size =1;
nx = 10;
ny = 10;
% time step size
dt = 1e-4;
% thermal diffusivity
alpha = 1.2;
% initial guess values for temperature
T = 300*ones(nx,ny);
% creating 2D mesh
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = 1/(nx-1);
dy = 1/(ny-1);
% tolerence for convergence of temperature

% initial conditions
for i = 1:nx
T(1,i) = 600;
T(ny,i) = 900;
end

for j = 1:ny
T(j,1) = 400;
T(j,nx) = 800;
end

T(1,1) = (600+400)/2;
T(1,nx) = (600+800)/2;
T(ny,1) = (400+900)/2;
T(ny,nx) = (900+800)/2;
%jaacobi method

k = 2*(dx^2+dy^2)/(dx^2*dy^2);
% convergence loop
tol = 1e-4;
T_initial = T;
T_old = T;
error = 100;
jacobi_iter = 1;

while(error > tol)
% space loop
for i = 2:nx-1
for j = 2:ny-1
% term_1 = (Told(i-1,j)+Told(i+1,j))/(dx^2);
% term_2 = (Told(i,j-1)+Told(i,j+1))/(dy^2);
% T(i,j) = (term_1+term_2)/k;
T(i,j) = (((T_old(i-1,j) + T_old(i+1,j))/(dx^2)) + ...
((T_old(i,j+1) + T_old(i,j-1))/(dy^2)))/k;
end
end
error = max(max(abs(T-T_old)));
T_old = T; % Updating Told values
jacobi_iter = jacobi_iter+1;
end

figure(1)
contourf(x,y,T)
figure(1)
colormap(jet)

A = contourf(x,y,T)
clabel(A)

set(gca,'YDir','reverse');

case_1 out put:-

gauss-siedel method octave code:-

% initial conditions
for i = 1:nx
T(1,i) = 600;
T(ny,i) = 900;
end

for j = 1:ny
T(j,1) = 400;
T(j,nx) = 800;
end

T(1,1) = (600+400)/2;
T(1,nx) = (600+800)/2;
T(ny,1) = (400+900)/2;
T(ny,nx) = (900+800)/2;
%gauss-siedel-method

k = 2*(dx^2+dy^2)/(dx^2*dy^2);
% convergence loop
tol = 1e-4;
T_initial = T;
T_old = T;
error = 100;
Gauss_iter = 1;

while(error > tol)
% space loop
for i = 2:nx-1
for j = 2:ny-1
% term_1 = (Told(i-1,j)+Told(i+1,j))/(dx^2);
% term_2 = (Told(i,j-1)+Told(i,j+1))/(dy^2);
% T(i,j) = (term_1+term_2)/k;
T(i,j) = (((T(i-1,j) + T_old(i+1,j))/(dx^2)) + ...
((T_old(i,j+1) + T(i,j-1))/(dy^2)))/k;
end
end
error = max(max(abs(T-T_old)));
T_old = T; % Updating Told values
Gauss_iter = Gauss_iter+1;
end

figure(2)
contourf(x,y,T)
figure(2)
colormap(jet)

A = contourf(x,y,T)
clabel(A)

set(gca,'YDir','reverse');

output of gauss -siedel-method:-

sor method case_3:-

code:-

% initial conditions
for i = 1:nx
T(1,i) = 600;
T(ny,i) = 900;
end

for j = 1:ny
T(j,1) = 400;
T(j,nx) = 800;
end

T(1,1) = (600+400)/2;
T(1,nx) = (600+800)/2;
T(ny,1) = (400+900)/2;
T(ny,nx) = (900+800)/2;

k = 2*(dx^2+dy^2)/(dx^2*dy^2);
% convergence loop
tol = 1e-4;
T_initial = T;
T_old = T;
error = 100;
sor_iter = 1;
omega = 1.5;

while(error > tol)
% space loop
for i = 2:nx-1
for j = 2:ny-1
term_1 = (T(i-1,j)+T_old(i+1,j))/(dx^2);
term_2 = (T(i,j-1)+T_old(i,j+1))/(dy^2);
% T(i,j) = (term_1+term_2)/k;
T(i,j) = ((1-omega)*T_old(i,j)) + (omega*(term_1+term_2))/k;
end
end
error = max(max(abs(T-T_old)));
T_old = T; % Updating Told values
sor_iter = sor_iter+1;
end

figure(3)
contourf(x,y,T)
figure(3)
colormap(jet)

A = contourf(x,y,T)
clabel(A)
set(gca,'YDir','reverse');

out-put:-

screen shorts of code:-

gauss-siedel-method:-

sor method:-

UNSTEADY STATE IMPLICITY FORM:-

CASE_1 [JACOBI METHOD]:-

clear all
close all
domain_size =1;
nx = 10;
ny = 10;
% time step size
dt = 0.01;
% thermal diffusivity
alpha = 0.4;
% initial guess values for temperature
T = 300*ones(nx,ny);
% creating 2D mesh
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = 1/(nx-1);
dy = 1/(ny-1);

% initial conditions
for i = 1:nx
T(1,i) = 600;
T(ny,i) = 900;
end

for j = 1:ny
T(j,1) = 400;
T(j,nx) = 800;
end

T(1,1) = (600+400)/2;
T(1,nx) = (600+800)/2;
T(ny,1) = (400+900)/2;
T(ny,nx) = (900+800)/2;
%jaacobi method

k1 = alpha*dt/(dx^2);
k2 = alpha*dt/(dy^2);
tol = 1e-4;
T_old = T;
Tprev = T_old;

jacobi_iter = 1;
nt = 1400;
for p = 1:nt
error = 9e9;

% convergence loop
while(error > tol)

% space loop
for i = 2:nx-1
for j = 2:ny-1
term_1 = (1+(2*k1)+(2*k2))^-1;
term_2 = k1*term_1;
term_3 = k2*term_1;
H = T_old(i-1,j)+T_old(i+1,j);
V = T_old(i,j-1)+T_old(i,j+1);
T(i,j) = (Tprev(i,j)*term_1)+(H*term_2)+(V*term_3);

end
end
error = max(max(abs(T_old-T)));
T_old = T; % Updating Told values
jacobi_iter = jacobi_iter+1;
end
Tprev = T;
end

figure(1)
colormap(jet)

A = contourf(x,y,T)
clabel(A)

set(gca,'YDir','reverse');

JACOBI OUTPUT:-

CASE_2 GAUSS-SIEDEL METHOD:-

% initial conditions
for i = 1:nx
T(1,i) = 600;
T(ny,i) = 900;
end

for j = 1:ny
T(j,1) = 400;
T(j,nx) = 800;
end

T(1,1) = (600+400)/2;
T(1,nx) = (600+800)/2;
T(ny,1) = (400+900)/2;
T(ny,nx) = (900+800)/2;
%Gauss siedel method

k1 = alpha*dt/(dx^2);
k2 = alpha*dt/(dy^2);
tol = 1e-4;
T_old = T;
Tprev = T_old;

Gauss_iter = 1;
nt = 1400;
for p = 1:nt
error = 9e9;

% convergence loop
while(error > tol)

% space loop
for i = 2:nx-1
for j = 2:ny-1
term_1 = (1+(2*k1)+(2*k2))^-1;
term_2 = k1*term_1;
term_3 = k2*term_1;
H = T(i-1,j)+T_old(i+1,j);
V = T(i,j-1)+T_old(i,j+1);
T(i,j) = (Tprev(i,j)*term_1)+(H*term_2)+(V*term_3);

end
end
error = max(max(abs(T_old-T)));
T_old = T; % Updating Told values
Gauss_iter = Gauss_iter+1;
end
Tprev = T;
end

figure(1)
colormap(jet)

A = contourf(x,y,T)
clabel(A)

set(gca,'YDir','reverse');

CODE SCREEN SHORTS:-

CASE-3 SOR METHOD:-

% initial conditions
for i = 1:nx
T(1,i) = 600;
T(ny,i) = 900;
end

for j = 1:ny
T(j,1) = 400;
T(j,nx) = 800;
end

T(1,1) = (600+400)/2;
T(1,nx) = (600+800)/2;
T(ny,1) = (400+900)/2;
T(ny,nx) = (900+800)/2;
%sor method

k1 = alpha*dt/(dx^2);
k2 = alpha*dt/(dy^2);
tol = 1e-4;
T_old = T;
Tprev = T_old;

sor_iter = 1;
nt = 1400;
omega = 1.1;
for p = 1:nt
error = 9e9;

% convergence loop
while(error > tol)

% space loop
for i = 2:nx-1
for j = 2:ny-1
term_1 = (1+(2*k1)+(2*k2))^-1;
term_2 = k1*term_1;
term_3 = k2*term_1;
H = T(i-1,j)+T_old(i+1,j);
V = T(i,j-1)+T_old(i,j+1);
T(i,j) = ((1-omega)*T_old(i,j))+omega*((term_1*Tprev(i,j))+(term_2*H)+(term_3*V));

end
end
error = max(max(abs(T_old-T)));
T_old = T; % Updating Told values
sor_iter = sor_iter+1;
end
Tprev = T;
end

figure(1)
colormap(jet)

A = contourf(x,y,T)
clabel(A)

set(gca,'YDir','reverse');

EXPLICITY FORM:-

clear all
close all
% input parameters
nx = 10; % number of nodes in x-direction
ny = 10; % number of nodes in y- direction
dt = 0.002; % time step
alpha = 1.4;
nt = 500;
% initial guess values for temperature
T = 300*ones(nx,ny);
% creating 2D mesh
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = 1/(nx-1);
dy = 1/(ny-1);
iter = 1;

% initial conditions
for i = 1:nx
T(1,i) = 600;
T(ny,i) = 900;
end

for j = 1:ny
T(j,1) = 400;
T(j,nx) = 800;
end

T(1,1) = (600+400)/2;
T(1,nx) = (600+800)/2;
T(ny,1) = (400+900)/2;
T(ny,nx) = (900+800)/2;
% creating a 2D mesh
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2)-x(1);
dy = y(2)-y(1);
T_old = T;
k1 = alpha*(dt/(dx^2));
k2 = alpha*(dt/(dy^2));
% time loop
for k = 1:nt
for i = 2:nx-1
for j = 2:ny-1
Term_1 = T_old(i,j);
Term_2 = (k1*(T_old(i+1,j))-(2*T_old(i,j))+T_old(i-1,j));
Term_3 = (k2*(T_old(i,j+1))-(2*T_old(i,j))+T_old(i,j-1));
T(i,j) = (Term_1+Term_2+Term_3);
end
end
end

T_old = T; % updating the Told values

figure(1)
colormap(jet)

A = contourf(x,y,T)
clabel(A)

set(gca,'YDir','reverse');

OUT PUT:-

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

Ship work like this
in your own course.

Talk to a career counsellor

Ship work like thisin your own course.

Talk to a career counsellor

Ship work like this
in your own course.