clc;

clear;

tic

format longG

rho=100;

cp=1000;

k=10;

T0=3;

TW=3;

TE=5;

s=10;

t=60000;

x=3;

dt=100;

dx=0.1;

sizet=t/dt;

sizex=x/dx;

Tmatrix=zeros(sizet,sizex);

for i=1

? ?Tmatrix(i,:)=T0;

end

%我想用TDMA解方程組，所以需要一堆東西去記錄ai1

%ap0始終不變，可以作為源項(xiàng)之一納入b向量

%ap1、ae1和aw1不隨時(shí)間變化但隨坐標(biāo)變化，所以應(yīng)該以向量形式給出：

%ap1=zeros(1,)

%隱格式算一算

for i=2:sizet

? ?%設(shè)置一個(gè)sizex*sizex的系數(shù)矩陣A，用于存放ai1

? ?A=zeros(sizex,sizex);

? ?%結(jié)果向量b

? ?b=zeros(sizex,1);

? ?for j=1:1

? ? ? ?aw1=2*k/dx;%均分網(wǎng)格，dxe=dxw=dx;

? ? ? ?ae1=k/dx;

? ? ? ?ap0=rho*cp*dx/dt;

? ? ? ?ap1=ap0+ae1+aw1;

? ? ? ?A(j,j)=ap1;

? ? ? ?A(j,j+1)=-ae1;

? ? ? ?b(j)=s*dx+ap0*Tmatrix(i-1,j)+aw1*TW;

? ?end

? ?for j=2:sizex-1

? ? ? ?ap0=rho*cp*dx/dt;

? ? ? ?aw1=k/dx;

? ? ? ?ae1=k/dx;

? ? ? ?ap1=ap0+aw1+ae1;

? ? ? ?A(j,j-1)=-aw1;

? ? ? ?A(j,j)=ap1;

? ? ? ?A(j,j+1)=-ae1;

? ? ? ?b(j)=s*dx+ap0*Tmatrix(i-1,j);

? ?end

? ?for j=sizex:sizex

? ? ? ?ap0=rho*cp*dx/dt;

? ? ? ?aw1=k/dx;

? ? ? ?ae1=2*k/dx;

? ? ? ?ap1=ap0+aw1+ae1;

? ? ? ?%ap1*Tmatrix(i,j)-aw1*Tmatrix(i,j-1)=ap0*Tmatrix(i-1,j)+s*dx+ae1*TE;

? ? ? ?A(j,j-1)=-aw1;

? ? ? ?A(j,j)=ap1;

? ? ? ?b(j)=s*dx+ap0*Tmatrix(i-1,j)+ae1*TE;

? ?end

? ?A;

? ?T1=A\b;

? ?Tmatrix(i,:)=T1';

end

? ?Ans=Tmatrix(sizet,:)'

? ?toc