1 function rechuandao() %
 2 
 3 Llist = ？？？
 4 N = 1000; % 空间点数
 5 M = 10000;
 6 alfa = ？？？ % 导热 / (密度*比热)
 7 lambda = 0.4; % 稳定条件
 8 %************
 9 zone=Llist(1)/N;%空间步长
10 z=0:zone:Llist(1);%空间点坐标数组
11 z=z';%矩阵转置
12 
13 n = lambda * (Llist(1)/N)^2 / alfa(1);
14 
15 TM = M * n;   % 总时间
16 
17 t=0:n:TM;    %时间点坐标数组
18 t=t';        %矩阵转置
19 
20 %计算初值和边值
21 T=zeros(N+1,M+1);  %赋温度矩阵初值
22 Ti=init_fun(z);    %各空间点处的初始条件
23 To=border_funo(t); %内边界条件
24 Te=border_fune(t); %外边界条件
25 T(:,1)=Ti;         %(初始条件 温度都为37)
26 T(1,:)=To;         %（边界条件x=0处温度为）
27 T(N+1,:)=Te;       %（边界条件x=L 处温度为）
28 
29 %用差分法求出温度T与杆长L、时间t的关系
30 for k=1:M  %时间点数
31     n=2;
32     while n<=N %空间点数
33     T(n,k+1)=lambda*(T(n+1,k)+T(n-1,k))+(-2*lambda+1)*T(n,k);
34     n=n+1;
35     end
36 end
37 %设置立体网格
38 for i=1:M+1
39     X(:,i)=z;
40 end
41 for j=1:N+1
42     Y(j,:)=t;
43 end
44 mesh(X,Y,T);
45 view([1 -1 1]);
46 xlabel('长度');
47 ylabel('时间');
48 zlabel('温度');
49 csvwrite('csv.csv',T)
50 
51 function y=init_fun(z)%初值条件
52 y=？？？
53 return
54 
55 function y=border_funo(t)%z=0的边界条件
56 y=？？？
57 return
58 
59 function y=border_fune(t)%z=L的边界条件
60 y= ？？？
61 return
62 %