This code performs heat transfer.

Tf=100; % Final Temp
Tint=1000; % Initial Temp
Tinf=20; % Water Temp
h=3000; % W/m^2*K
rho=5550; %kg/m^3
dx=0.005; %m
L=.25/2;%m
num=(L/dx);
T=zeros(num,1);
T(:,1)=Tint;
Cp=600;
k=27.5;
Alpha=(k/(rho*Cp));
Bi=((h*dx)/k);
dt=(dx^2)/(2*Alpha*(1+Bi));
Fo=((Alpha*dt)/(dx^2));
i=1;
while T(num,i)>Tf
for j=1:num
if j==1
T(j,i+1)=(2*Fo(j)*Bi(j))*Tinf+(2*Fo(j))*T(j+1,i)+(1-(2*Fo(j)*Bi(j))-(2*Fo(j)))*T(j,i);
elseif j==num
T(j,i+1)=2*Fo(j)*T((j-1),i)+(1-(2*Fo(j)))*T(j,i);
else
T(j,i+1)=(Fo(j))*(T((j-1),i)+T((j+1),i))+(1-(2*Fo(j)))*T(j,i);
end
if (20<=(T(j,i)))&&((T(j,i))<600)
Cp(j)=425+(.773)*(T(j,i))-.00169*((T(j,i))^2)+.00000222*((T(j,i))^3);
elseif (600<=(T(j,i)))&&(T(j,i))<(735)
Cp(j)=666+(13002/(738-(T(j,i))));
elseif (735<=(T(j,i)))&&(T(j,i))<(900)
Cp(j)=545+(17820/((T(j,i))-731));
elseif (900<=(T(j,i)))&&(T(j,i))<=(1000)
Cp(j)=600;
end
if(20<=(T(j,i)))&&((T(j,i))<800)
k(j)=54-.0333*(T(j,i));
elseif (800<=(T(j,i)))&&(T(j,i))<=(1000)
k(j)=27.3;
end
Alpha=(k(j)/(rho*Cp(j)));
Bi(j)=((h*dx)/k(j));
Fo(j)=((Alpha*dt)/(dx^2));
end
i=i+1;
end
Index exceeds the number of array elements. Index must not exceed 1.

댓글 수: 3

Torsten
Torsten 2023년 4월 25일
편집: Torsten 2023년 4월 25일
Sorry, but your code is full of errors. Maybe it makes more sense if you explain what you are trying to do.
M Miller
M Miller 2023년 4월 26일
This is part of the code to find how long it takes a slab of steele to cool from 1200° C to 200° C by using water jets at 40° C
Torsten
Torsten 2023년 4월 26일
Seems you try to solve the transient heat conduction equation in 1d.
MATLAB's "pdepe" is the correct tool to do this.

댓글을 달려면 로그인하십시오.

답변 (1개)

Cris LaPierre
Cris LaPierre 2023년 4월 25일

0 개 추천

Your code is written expecting Fo and Bi to be a vectors, but they are scalars, so only ever has 1 element. You get this error any time you try to index using a value greater than 1.
The easiest solution based on the information we have is to not index these 2 variables.
Tf=200; % Final Temp
Tint=1200; % Initial Temp
Tinf=40; % Water Temp
h=4000; % W/m^2*K
rho=7850; %kg/m^3
dx=0.005; %m
L=.5/2;%m
num=(L/dx)+1;
T=zeros(num,1);
T(:,1)=Tint;
Cp=650; % initial Cp
k=27.3; % Initial k
% used to solve for dt
Alpha=(k/(rho*Cp));
Bi=((h*dx)/k);
dt=(dx^2)/(2*Alpha*(1+Bi));
Fo=((Alpha*dt)/(dx^2));
i=1;
while T(num,i)>Tf
for j=1:num
if j==1
T(j,i+1)=(2*Fo*Bi)*Tinf+(2*Fo)*T(j+1,i)+(1-(2*Fo*Bi)-(2*Fo))*T(j,i);
elseif j==num
T(j,i+1)=2*Fo*T((j-1),i)+(1-(2*Fo))*T(j,i);
else
T(j,i+1)=(Fo)*(T((j-1),i)+T((j+1),i))+(1-(2*Fo))*T(j,i);
end
% Cp and k need to be defined each time it finds a temperature since they are both temperature dependent
if (20<=(T(j,i)))&&((T(j,i))<600)
Cp(j)=425+(.773)*(T(j,i))-.00169*((T(j,i))^2)+.00000222*((T(j,i))^3);
elseif (600<=(T(j,i)))&&(T(j,i))<(735)
Cp(j)=666+(13002/(738-(T(j,i))));
elseif (735<=(T(j,i)))&&(T(j,i))<(900)
Cp(j)=545+(17820/((T(j,i))-731));
elseif (900<=(T(j,i)))&&(T(j,i))<=(1200)
Cp(j)=650;
end
if(20<=(T(j,i)))&&((T(j,i))<800)
k(j)=54-.0333*(T(j,i));
elseif (800<=(T(j,i)))&&(T(j,i))<=(1200)
k(j)=27.3;
end
%These values are recalculated each time as well based on the Cp and k
Alpha=(k(j)/(rho*Cp(j)));
Bi=((h*dx)/k(j));
Fo=((Alpha*dt)/(dx^2));
end
i=i+1;
end
ttot=(i-1)*dt
ttot = 826.5942

카테고리

도움말 센터File Exchange에서 Thermal Analysis에 대해 자세히 알아보기

태그

질문:

2023년 4월 25일

편집:

2023년 4월 28일

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by