2D Heat conduction steady and transient

steady state:

Since dx and dy are same, the equation reduces to a simple form.

clear all
close all
clc



nx=11;
ny=11;
x=linspace(0,1,nx);
y=linspace(0,1,ny);

method=\'Enter which method to use:\';
method = input(method);



c=zeros(1,nx);
t = meshgrid(c);


    t(1:nx,1) = 400; %T_left
    t(1,1:ny) = 600; %T_top
    t(nx,1:ny) = 900; %T_bottom
    t(1:ny,ny) = 800; %T_right

%Corner temperatures
t(1,1) = 500;
t(nx,1) = 650;
t(nx,ny) = 850;
t(1,ny) = 700;



told=t;
tinit=t;
tol = 1e-4;
error = 400;

jacob = 1;
gauss = 1;


if( method ==1)
while(error > tol)
    for i=2:nx-1
    for j=2:ny-1
t(i,j) = 0.25 * (told(i-1,j) + told(i+1,j) + told(i,j-1) + told(i,j+1));
    end
    end

error = max(max(abs(t-told)));
told = t;

contour(x,y,t);
xlabel(\'x\');
ylabel(\'y\');
colorbar;
heading=sprintf(\'Iteration number = %d\',jacob);
title(heading);
pause(0.03);
jacob = jacob+1;


end

end

if( method ==2)
while(error > tol)
    for i=2:nx-1
    for j=2:ny-1
t(i,j) = 0.25 * (t(i-1,j) + told(i+1,j) + t(i,j-1) + told(i,j+1));
    end
    end

error = max(max(abs(t-told)));
told = t;


contour(x,y,t);
xlabel(\'x\');
ylabel(\'y\');
colorbar;
heading=sprintf(\'Iteration number = %d\',gauss);
title(heading);
pause(0.03);
gauss = gauss+1;

end
end

steady stat jacobi steady state gauss

Explanation:

The gauss seidel method converges faster than the jacobi method, because the gauss method uses updated values for every iteration.

unsteady explicit:

The explicit function uses old values which need to be updated for every iteration. This is in a way similar to the steady state response.

clear all
close all
clc

nx=11;
ny=11;
x=linspace(0,1,nx);
dx=1/(nx-1);
y=linspace(0,1,ny);
dy=1/(ny-1);
alpha = 0.01;
dt=0.25;

method=\'Enter which method to use:\';
method = input(method);



c=zeros(1,nx);
t = meshgrid(c);


    t(1:nx,1) = 400; %T_left
    t(1,1:ny) = 600; %T_top
    t(nx,1:ny) = 900; %T_bottom
    t(1:ny,ny) = 800; %T_right

%Corner temperatures
t(1,1) = 500;
t(nx,1) = 650;
t(nx,ny) = 850;
t(1,ny) = 700;


told=t;
tprev=t;
tinit=t;
tol = 1e-4;
error = 400;
errort=400;


jacob = 1;
gauss = 1;

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

if(method==1)
while(errort>tol)
       
        for i=2:nx-1
            for j=2:ny-1
        t(i,j) = (tprev(i,j)/den) + ( (k1/den)*(told(i,j-1)+told(i,j+1)) ) + ( (k2/den)*(told(i+1,j)+told(i-1,j)) );

            end
        end
        
       

        errort=max(max(abs(t-tprev)));

        tprev= t;
        told = t;
        contour(x,y,t);
        
        xlabel(\'x\');
        ylabel(\'y\');
        heading=sprintf(\'Iteration number = %d\',jacob);
        title(heading);
        colorbar;
drawnow;

       
        jacob = jacob+1;
end
end




if(method==2)
while(errort>tol)
       
        for i=2:nx-1
            for j=2:ny-1
        t(i,j) = (tprev(i,j)/den) + ( (k1/den)*(t(i,j-1)+told(i,j+1)) ) + ( (k2/den)*(told(i+1,j)+t(i-1,j)) );

            end
        end
        
       

        errort=max(max(abs(t-tprev)));

        tprev= t;
        told = t;
        contour(x,y,t);
        
        xlabel(\'x\');
        ylabel(\'y\');
        heading=sprintf(\'Iteration number = %d\',gauss);
        title(heading);
        colorbar;
drawnow;

       
        gauss = gauss+1;
end
end

unsteady explicit jacobi unsteady explicit gauss

Explanation:

The CFL value is set in a way such that the value converges. While updating values inside the space loop, it is found that the value converged within one time interval, so values are being updated inside the time loop, and the contour was plotted.

unsteady implicit:

Present values are found using the predicted values of the next iteration, and then values are updated.

clear all
close all
clc

nx=11;
ny=11;
x=linspace(0,1,nx);
dx=1/(nx-1);
y=linspace(0,1,ny);
dy=1/(ny-1);
alpha = 0.01;
dt=0.25;

method=\'Enter which method to use:\';
method = input(method);



c=zeros(1,nx);
t = meshgrid(c);


    t(1:nx,1) = 400; %T_left
    t(1,1:ny) = 600; %T_top
    t(nx,1:ny) = 900; %T_bottom
    t(1:ny,ny) = 800; %T_right

%Corner temperatures
t(1,1) = 500;
t(nx,1) = 650;
t(nx,ny) = 850;
t(1,ny) = 700;


told=t;
% tprev=t;
tinit=t;
tp1 = t;
tol = 1e-4;
error = 400;
errort=400;


jacob = 1;
gauss = 1;

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

if(method==1)
while(errort>tol)
        for i=2:nx-1
            for j=2:ny-1
        tp1(i,j) = (t(i,j)/den) + ( (k1/den)*(told(i,j-1)+told(i,j+1)) ) + ( (k2/den)*(told(i+1,j)+told(i-1,j)) );
                                told=tp1;
            end
        end

        errort=max(max(abs(tp1-t)));

        t= tp1;
%         told = tp1;
%         told = t;
        contour(x,y,t);
        
        xlabel(\'x\');
        ylabel(\'y\');
        heading=sprintf(\'Iteration number = %d\',jacob);
        title(heading);
        colorbar;
drawnow;

       
        jacob = jacob+1;
end
end


if(method==2)
while(errort>tol)
        for i=2:nx-1
            for j=2:ny-1
        tp1(i,j) = (t(i,j)/den) + ( (k1/den)*(tp1(i,j-1)+told(i,j+1)) ) + ( (k2/den)*(told(i+1,j)+tp1(i-1,j)) );
                                told=tp1;
            end
        end

        errort=max(max(abs(tp1-t)));

        t= tp1;
%         told = tp1;
%         told = t;
        contour(x,y,t);
        
        xlabel(\'x\');
        ylabel(\'y\');
        heading=sprintf(\'Iteration number = %d\',gauss);
        title(heading);
        colorbar;
drawnow;

       
        gauss = gauss+1;
end
end

unsteady implicit jacobi unsteady implicit gauss

Explanation:

Based on the diffusivity values and time intervals, the time taken to attain steady state condition varies. For optimum case, I have taken these values such that CFL can be near to 0.25,