一维非稳态导热微分方程的数值求解

%%详细推导求解过程请移步 微信公众号 CFD入门指南 %%

clc
clear
L=0.2;     %杆长
t=10;      %计算时长
T1=20;     %杆初始温度
T2=80;     %环境温度
rou=8900;   %密度
c=390;    %比热容
lan=377;    %导热系数
h=10000;   %与环境的换热系数
dx=0.001;   %空间步长
dt=0.0001;  %时间步长
m=t/dt;    %时间单元数量
n=L/dx;    %空间单元数量
a=lan/rou/c;    
Fo=a*dt/dx/dx;   %傅里叶数
Bi=h*dx/lan;     %毕渥数

v1=ones（m+1n+1）*T1;     %初始时杆上均匀温度为T1

if Fo>0.5 && （1-2*Fo*Bi-2*Fo<0） 
    fprintf（‘无法收敛，调整步长‘）;    %判断是否满足收敛条件
else
    for i=1:m
        for j=2:n
             v1（i+1j）=