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

Steady state analysis Transient State Analysis

MMalluvajhula Manogna13 Feb 2019
MATLAB

2D heat conduction equation can be given as

1.Transient Equation

dT/dt + alpha((d^2T/dx^2)+(d^2T/dy^2))=0

2. Steady State Equaton

(d^2T/dx^2)+(d^2T/dy^2) = 0

Solutions to above equations can be given in both implicit and explicit forms for transient state and steady state.

1.Transient state

a. Implicit Method

T(i,j) = Tprev(i,j) + (k1*(Told(i-1,j) + Told(i+1,j))+(k2*(Told(i,j-1)+Told(i,j+1))

b. Explicit Method

T(i,j) = Told(i,j) + (k1*(Told(i-1,j) + Told(i+1,j)-(2*Told(i,j))+(k2*(Told(i,j-1)+Told(i,j+1)-(2*Told(i,j))

where k1 = (alpha*dt)/dx^2

k2 = (alpha*dt)/dy^2

2. Steady State

T(i,j)=((Told(i-1,j)+Told(i+1,j))+Told(i,j-1)+Told(i,(j+1)))*0.25

i. Jacobi Method

In jacobi method we assume old values of T

ii. Gauss Method

In gauss method we assume updated values at left and bottom points

iii. SOR Method

In sor method we consider the constant f

where f =2/(1+sin(pi*dx))

and equatiion is T(i,j)=((1-f)*Told(i,j))+((Told(i,j)+(k1*(Told(i+1,j)-(2*Told(i,j))+T(i-1,j)))+(k2*(Told(i,j+1)-(2*Told(i,j))+T(i,j-1))))*f)

A code is written in MATLAB to solve all the above iterative methods

Steady State Equation

%Steady State Equation
clear all
close all
clc

nx=5;
ny=5;
T = 300*ones(nx,ny);

T(nx,:) = 900;
T(1,:) = 600;
T(:,1) = 400;
T(:,ny) = 800;
Told = T;
x=linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2)-x(1);
dy = y(2)-y(1);

%1. jacobi_method
%2.gauss_iter
%3.SOR_Method

iter_solver = input(\'Enter the iter_solver value: \')
tol = 1e-4;
[X,Y] = meshgrid(x,y);
error = 9e9;
%jacobi_method
tic
if(iter_solver==1)
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))+Told(i,j-1)+Told(i,(j+1)))*0.25;
   end
   end
   error = max(abs(Told-T));
   Told = T;
      figure(1)
   [C,h] = contourf(X,Y,T);
   clabel(C,h)
   xlabel(\'x\')
   ylabel(\'y\')  
   jacobi_iter = jacobi_iter+1;
   end
   toc
  end
  
  %%Gauss Siedel Method
  tic
if(iter_solver==2)
gauss_iter =1;
while(error>tol)
for i=2:(nx-1)
  for j=2:(ny-1)
    T(i,j)=((T(i-1,j)+Told(i+1,j))+T(i,j-1)+Told(i,(j+1)))*0.25;
   end
   end
   error = max(abs(Told-T));
   Told = T;
   figure(2)
   [C,h] = contourf(X,Y,T);
   clabel(C,h)
   xlabel(\'x\')
   ylabel(\'y\')  
   gauss_iter = gauss_iter+1;
   end
   toc
  end
  
  
  %SOR_Method
    tic
if(iter_solver==3)
sor_iter =1;
a = 2/(1+sin(pi*dx));
while(error>tol)
for i=2:(nx-1)
  for j=2:(ny-1)
    T(i,j)=((1-a)*(Told(i,j)))+(((T(i-1,j)+Told(i+1,j))+T(i,j-1)+Told(i,(j+1)))*0.25*a);
   end
   end
   error = max(abs(Told-T));
   Told = T;
   figure(3)
   [C,h] = contourf(X,Y,T);
   clabel(C,h)
    xlabel(\'x\')
   ylabel(\'y\')  
   sor_iter = sor_iter+1;
   end
   toc
  end

Graphs for the above equations are

1.Jacobi Method

output:

2. Gauss Siedel Method

output:

3. SOR Method

output

Transient state equation

The main code consists of functions of different iterative methods

clear all
close all
clc

for iter_solver=1:6
 if (iter_solver==1)
 T = jacobi_1(iter_solver);
 elseif (iter_solver==2)
 T = gauss_1(iter_solver);
 elseif (iter_solver==3)
 T = sor_1(iter_solver);
 elseif (iter_solver==4)
 T = jacobi_2(iter_solver);
 elseif (iter_solver==5)
 T = gauss_2(iter_solver);
 elseif (iter_solver==6)
 T = sor_2(iter_solver);
 end
 end

i. Explicit jacobi function

function result1 =  jacobi_1(iter_solver)
% material = coper

nx = 5;
ny = 5;
T = 300*ones(nx,ny);

T(nx,:) = 900;
T(1,:) = 600;
T(:,1) = 400;
T(:,ny) = 800;
Told = T;
Tprev = T;
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2)-x(1);
dy = y(2)-y(1);
error = 1;
tol = 1e-4;
[X,Y] = meshgrid(x,y);

a= 1.11*10^-4;
dt = 0.01;
k1 = (a*dt)/(dx*dx);
k2 = (a*dt)/(dy*dy);
tic
jacobi_iter =1;
while (error>tol)
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
   error = max(abs(Told-T));
   Told = T;
      jacobi_iter = jacobi_iter+1;
      figure(1)
   [C,h] = contourf(X,Y,T);
   clabel(C,h)
   xlabel(\'x\')
   ylabel(\'y\')
   
   end
   toc
   simulationtime=toc;
   title(sprintf(\'Explicit Jacobi method and Simulation time = %f \', simulationtime));
   
  end

Graph

ii. Explicit Gauss siedel method

function result2 =  gauss_1(iter_solver)
% material = coper

nx = 5; 
ny = 5;
T = 300*ones(nx,ny);

T(nx,:) = 900;
T(1,:) = 600;
T(:,1) = 400;
T(:,ny) = 800;
Told = T;
Tprev = T;
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2)-x(1);
dy = y(2)-y(1);
error = 9e9;
tol = 1e-4;
[X,Y] = meshgrid(x,y);

a= 1.11*(10^-4);
dt = 0.01;
k1 = (a*dt)/(dx*dx);
k2 = (a*dt)/(dy*dy);
tic
gauss_iter =1;
while (error>tol) 

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))+T(i-1,j)))+(k2*(Told(i,j+1)-(2*Told(i,j))+T(i,j-1)));
   end
   end
   error = max(abs(Told-T));
   Told = T;
    gauss_iter = gauss_iter+1;
      figure(2)
   [C,h] = contourf(X,Y,T);
   clabel(C,h)
   xlabel(\'x\')
   ylabel(\'y\')
   pause(0.03)
  
   end 
   toc
   simulationtime=toc;
   title(sprintf(\'Explicit Gauss method and Simulation time = %f \', simulationtime));
   
  end

Graph

iii. Explicit SOR method

function result3 =  sor_1(iter_solver)
% material = coper


nx = 5;
ny = 5;
T = 300*ones(nx,ny);

T(nx,:) = 900;
T(1,:) = 600;
T(:,1) = 400;
T(:,ny) = 800;
Told = T;
Tprev = T;
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2)-x(1);
dy = y(2)-y(1);
error = 9e9;
tol = 1e-4;
[X,Y] = meshgrid(x,y);

a= 1.11*10^-4;
f =2/(1+sin(pi*dx));
dt = 0.01;
k1 = (a*dt)/(dx*dx);
k2 = (a*dt)/(dy*dy);
tic
sor_iter =1;
while (error>tol)
for i=2:(nx-1)
  for j=2:(ny-1)
    T(i,j)=((1-f)*Told(i,j))+((Told(i,j)+(k1*(Told(i+1,j)-(2*Told(i,j))+T(i-1,j)))+(k2*(Told(i,j+1)-(2*Told(i,j))+T(i,j-1))))*f);
   end
   end
   error = max(abs(Told-T));
   Told = T;
    sor_iter = sor_iter+1;
      figure(3)
   [C,h] = contourf(X,Y,T);
   clabel(C,h)
   xlabel(\'x\')
   ylabel(\'y\')  
   end
   toc
   simulationtime=toc;
   title(sprintf(\'Explicit SOR method and Simulation time = %f \', simulationtime));
  
  end

Graph

iv. Implicit jacobi method

function result4 =  jacobi_2(iter_solver)
% material = coper


nx = 5;
ny = 5;
T = 300*ones(nx,ny);

T(nx,:) = 900;
T(1,:) = 600;
T(:,1) = 400;
T(:,ny) = 800;
Told = T;
Tprev = T;
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2)-x(1);
dy = y(2)-y(1);
error = 9e9;
tol = 1e-4;
[X,Y] = meshgrid(x,y);

a= 1.11*10^-4;
dt = 0.01;
k1 = (a*dt)/(dx*dx);
k2 = (a*dt)/(dy*dy);
b = 1/(1+(2*k1)+(2*k2));
c = b*k1;
d= b*k2;
tic
jacobi_iter =1;
while (error>tol)
for i=2:(nx-1)
  for j=2:(ny-1)
    T(i,j)=(b*Tprev(i,j))+(c*(Told(i+1,j)-(2*Told(i,j))+Told(i-1,j)))+(d*(Told(i,j+1)-(2*Told(i,j))+Told(i,j-1)));
   end
   end
   error = max(abs(Told-T));
   Told = T;
    jacobi_iter = jacobi_iter+1;
      figure(4)
   [C,h] = contourf(X,Y,T);
   clabel(C,h)
   xlabel(\'x\')
   ylabel(\'y\')
   
   end
   toc
   simulationtime=toc;
   title(sprintf(\'Implicit Jacobi method and Simulation time = %f \', simulationtime));
  end

Graph

v. Implicit Gauss Siedel method

function result5 =  gauss_2(iter_solver)
% material = coper


nx = 5;
ny = 5;
T = 300*ones(nx,ny);

T(nx,:) = 900;
T(1,:) = 600;
T(:,1) = 400;
T(:,ny) = 800;
Told = T;
Tprev = T;
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2)-x(1);
dy = y(2)-y(1);
error = 9e9;
tol = 1e-4;
[X,Y] = meshgrid(x,y);

a= 1.11*10^-4;
dt = 0.01;
k1 = (a*dt)/(dx*dx);
k2 = (a*dt)/(dy*dy);
b = 1/(1+(2*k1)+(2*k2));
c = b*k1;
d= b*k2;
tic
gauss_iter =1;
while (error>tol)
for i=2:(nx-1)
  for j=2:(ny-1)
    T(i,j)=(b*Tprev(i,j))+(c*(Told(i+1,j)-(2*Told(i,j))+T(i-1,j)))+(d*(Told(i,j+1)-(2*Told(i,j))+T(i,j-1)));
   end
   end
   error = max(abs(Told-T));
   Told = T;
  gauss_iter = gauss_iter+1; 
      figure(5)
   [C,h] = contourf(X,Y,T);
   clabel(C,h)
   xlabel(\'x\')
   ylabel(\'y\')
   
   end
   toc
   simulationtime=toc;
   title(sprintf(\'Implicit Gauss method and Simulation time = %f \', simulationtime));
  end

Graph

vi. Implicit SOR Method

function result6 =  sor_2(iter_solver)
% material = coper


nx = 5;
ny = 5;
T = 300*ones(nx,ny);

T(nx,:) = 900;
T(1,:) = 600;
T(:,1) = 400;
T(:,ny) = 800;
Told = T;
Tprev = T;
x = linspace(0,1,nx);
y = linspace(0,1,ny);
dx = x(2)-x(1);
dy = y(2)-y(1);
error = 9e9;
tol = 1e-4;
[X,Y] = meshgrid(x,y);

a= 1.11*10^-4;
f =2/(1+sin(pi*dx));
dt = 0.01;
k1 = (a*dt)/(dx*dx);
k2 = (a*dt)/(dy*dy);
b = 1/(1+(2*k1)+(2*k2));
c = b*k1;
d= b*k2;
tic
sor_iter =1;
while (error>tol)
for i=2:(nx-1)
  for j=2:(ny-1)
    T(i,j)=((1-f)*Told(i,j))+(((b*Tprev(i,j))+(c*(Told(i+1,j)-(2*Told(i,j))+T(i-1,j)))+(d*(Told(i,j+1)-(2*Told(i,j))+T(i,j-1))))*f);
   end
   end
   error = max(abs(Told-T));
   Told = T;
    sor_iter = sor_iter+1;
      figure(6)
   [C,h] = contourf(X,Y,T);
   clabel(C,h)
   xlabel(\'x\')
   ylabel(\'y\')
   
   end
   toc
   simulationtime=toc;
   title(sprintf(\'Implicit SOR method and Simulation time = %f \', simulationtime));
  end

Graph

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.