1

First of all I am just a starter in Matlab, so the program is not optimized. For now I just want to get it running. This is a 3D heat problem with boundary conditions. Fixed temperature from the north, insulated from the bottom, and we have convection on the 4 surrounding sides with the same h. when I run it it does not give me the right convergences. I have the analytical and the solidworks approximations, so I know the results I get in Matlab are wrong. Could someone advise me please. Your help is appreciated.

==============================

%===Solution of the problem of 2D in the bookmarks:Case2a=====
clc;
close all;
clear all;

DX=0.01; % step size
DY=0.01;
DZ=0.01;
Lx= 0.6; %length along x-axis in cmn
Ly=0.4;   %length along y-axis
Lz=0.3;   %length along z-axis
Nx=Lx/DX+1;    %Number of nodes along x-axis
Nx=floor(Nx)+1; % I shouldn't add 1. I should get it automatically 
Ny=(Ly/DY+1);    %Number of nodes along y-axis
Ny=floor(Ny);
Nz=(Lz/DZ+1);
Nz=floor(Nz)+1;

K=2;
h=500;
T_inf=20;

X=0:DX:Lx;
Y=0:DY:Ly;
Z=0:DZ:Lz;

T=zeros(Nx,Ny,Nz);                                                            

% I added imaginary nodes outside then it gave me numbers

tol=1;
s=0; %s=2000 could be set as the maximum number of allowed iteration for 
example
err=1; 
T_old=1;

while err >=tol %&& s<2001
    s=s+1;

%------------boundary conditions------------------------------------------

%====North temperature (Boundary)
T(:,Ny,:)=90;  
%===South Boundary "insulation" ==============
for k=1:Nz
    for i=2:Nx-1
   T(i,1,k)=0.25*[2*T(i,2,k)+T(i-1,1,k)+T(i+1,1,k)];
  %T(i-1,1,k)=4* T(i,1,k)-2*T(i,2,k)-T(i+1,1,k);
    end
end
%================Back Boundary "Convection" ================
for j=2:Ny-1 %changed
    for i=2:Nx-1
  T(i,j,1)=0.5*K/(h*DX+2*K)*[2*T(i,j,2)+T(i-1,j,1)+T(i+1,j,1)+2*h*DX*T_inf/K];
    end
end
%================Front Boundary "Convection" ================
for j=2:Ny-1
    for i=2:Nx-1
  T(i,j,Nz)=0.5*K/(h*DX+2*K)*[2*T(i,j,Nz-1)+T(i-1,j,Nz-1)+T(i+1,j,Nz-1)+2*h*DX*T_inf/K];
    end
end

%================West Boundary "Convection" ================
for j=2:Ny-1
for k=2:Nz-1
  T(1,j,k)=0.5*K/(h*DX+2*K)*[2*T(2,j,k)+T(1,j,k-1)+T(1,j,k+1)+2*h*DX*T_inf/K];
    end
end
%================East Boundary "Convection" ================
for j=2:Ny-1
    for k=2:Nz-1
  T(Nx,j,k)=0.5*K/(h*DX+2*K)*[2*T(Nx-1,j,k)+T(Nx,j,k-1)+T(Nx,j,k+1)+2*h*DX*T_inf/K];
    end
end



% The interior nodes
for k = 2:Nz-1 % I think it needs to go loop by loop(pairs of loops)
    for j = 2:Ny-1
        for i=2:Nx-1
    T(i,j,k)=0.25*(T(i-1,j)+T(i+1,j)+T(i,j-1)+T(i,j+1));
        end
    end
end
err=max(max(abs(T-T_old)));
T_old=T;
end

T=rot90(T)
  • maybe more a stackoverflow.com question, but be aware "please debug my code, it doesn't work" is not a good way to get help there either. – agentp Nov 28 '17 at 03:09
  • Actually this question was asked there first, and they were the ones who referred me to this section. – Rachid Brah Nov 28 '17 at 13:16
  • well, you had only one comment to that effect. https://stackoverflow.com/q/47512200/1004168 Its really not clear where it goes since you have not provided any explanation of the problem you are having. – agentp Nov 28 '17 at 13:27
  • the 3D matrix I get when I run the code above gives me 0 at the four vertical edges of the cube. That is my issue – Rachid Brah Nov 28 '17 at 15:33
  • I think https://compsci.stackexchange.com may be a better match for your question. – peterh Nov 29 '17 at 18:13

0 Answers0