二维非稳态导热MATLAB程序

clc;
clear;
options={'长宽L','空间点数N','扩散系数alfa','稳定条件的值lambda(取值必须小于0.5)','时间','时间步长'};
topic='seting';
lines=1;
def={'0.1','10','0.05','0.5','10','0.001'};
h=inputdlg(options,topic,lines,def);
L=eval(h{1});
N=eval(h{2});
alfa=eval(h{3});
lambda=eval(h{4});%lambda的值必须小于0.5
t=eval(h{5});
timestep=eval(h{6});
%***************************************************
h=L/N;%空间步长
x=0:h:L;
x=x';
y=0:h:L;
y=y';
tao=(alfa*timestep)/(h*h);
T=zeros(N+1,N+1,t);
% 初始条件
T(:,:,1)=20;

for i=1:1:0.1/0.001;

for k=2:N
m=2;

while m<=N
% 外部脚点

T(N+1,1,i+1)=4*tao*h*(T(N+1,2,i)+T(N,1,i))+(1-4*tao)*T(N+1,1,i)+0.1;
T(1,N+1,i+1)=4*tao*h*(T(1,N,i)+T(2,N+1,i))+(1-4*tao)*T(1,N+1,i)+0.1;
T(N+1,N+1,i+1)=4*tao*h*(T(N,N+1,i)+T(N+1,N,i))+(1-4*tao)*T(N+1,N+1,i)+0.1;
T(1,1,i+1)=4*tao*h*(T(1,2,i)+T(2,1,i))+(1-4*tao)*T(1,1,i)+0.1;
% 平直边界
T(1,k,i+1)=tao*h*(2*T(2,k,i)+T(1,k-1,i)+T(1,k+1,i))+(1-4*tao)*T(1,k,i)+0.2;
T(N+1,k,i+1)=tao*h*(2*T(N,k,i)+T(N+1,k-1,i)+T(N+1,k+1,i))+(1-4*tao)*T(N+1,k,i)+0.2;
T(m,N+1,i+1)=tao*h*(2*T(m,N,i)+T(m-1,N+1,i)+T(m+1,N+1,i))+(1-4*tao)*T(m,N+1,i)+0.2;
T(m,1,i+1)=tao*h*(2*T(m,2,i)+T(m-1,2,i)+T(m+1,2,i))+(1-4*tao)*T(m,1,i)+0.2;

T(m,k,i+1)=tao*(T(m+1,k,i)+T(m-1,k,i)+T(m,k+1,i)+T(m,k-1,i))+(1-4*tao)*T(m,k,i)+0.4;
m=m+1;

end;


end
end