哈尔滨工程大学传热学大作业数值计算matlab程序内容

传热学作业数值计算
1 
程序内容： 数值计算matlab程序内容：>> tw1=10; % 赋初值
赋初值 tw2=20; 
c=1.5; 
p2=20; 
p1=c*p2; 
L2=40; 
L1=c*L2; 
deltaX=L2/p2; 
a=p2+1; 
b=p1+1; 
=ones(a,b)*5; 
m1=ones(a,b); 
m1(a,2:b-1)=zeros(1,b-2); 
m1(2:a,1)=zeros(a-1,1); 
m1(2:a,b)=zeros(a-1,1); 
m1(1,:)=ones(1,b)*2; 
k=0; 
max1=1.0; 
tn= ; 
while(max1>1e-6) 
max1=0; 
k=k+1; 
for i=1:1:a 
for j=1:1:b 
m=m1(i,j); 
n= (i,j); 
switch m 
case 0 
tn(i,j)=tw1; 
case 1 
tn(i,j)=0.25*(tn(i,j+1)+tn(i,j-1)+tn(i+1,j)+tn(i-1,j)); 
case 2 
tn(i,j)=tw1+tw2*sin(pi*(j-1)/(b-1)); 
end 
er=abs(tn(i,j)-n); 
if er>max1 
max1=er; 
end 
end 
end 
=tn; 
end 
k 
max1 
t2=ones(a,b); %求解析温度场
求解析温度场
for i=a:-1:1 
for j=1:1:b 
y=deltaX*(a-i); 
x=deltaX*(j-1); 
t2(i,j)=tw1+tw2*sin(pi*x/L1)*(sinh(pi*y/L1))/(sinh(pi*L2/L1)); 
end 
end 
t2 
迭代次数k =706 
数值解温度场
数值解每次迭代的最大误差max1 =9.8531e-07 
解析温度场 t2 
解析温度场
取第11行的解析解和数值解的点
行的解析解和数值解的点
行的解析解的直线，散点为其数值解的点 曲线为第11行的解析解的直线，散点为其数值解的点
解析解 第11行的误差=[数值解(11行) –解析解(11行)]/解析解
数值温度场图像数值温度场图像 解析温度场图像解析温度场图像
数值解与解析解的误差数值解与解析解的误差
数值计算matlab 程序内容：程序内容： >> tw1=10; tw2=20; c=1.5; p2=20; p1=c*p2; L2=20; 
deltaX=L2/p2; L1=c*L2; a=p2+1; b=p1+1; 
=ones(a,b)*5; 
m1=ones(a,b); 
m1(a,2:b-1)=zeros(1,b-2); 
m1(2:a,1)=zeros(a-1,1); 
m1(2:a,b)=zeros(a-1,1); 
m1(1,:)=ones(1,b)*2; 
k=0; 
max1=1.0; 
tn= ; 
while(max1>1e-6) 
max1=0; 
k=k+1; 
for i=1:1:a 
for j=1:1:b 
m=m1(i,j); 
n= (i,j); 
switch m 
case 0 
tn(i,j)=tw1; 
case 1 
tn(i,j)=0.25*(tn(i,j+1)+tn(i,j-1)+tn(i+1,j)+tn(i-1,j)); 
case 2 
tn(i,j)=tw2; 
end 
er=abs(tn(i,j)-n); 
if er>max1 
max1=er; 
end 
end 
end 
=tn; 
end 
k 
max1 
tx=ones(a,b); 
for i=1:1:a 
for j=1:1:b 
y=(a-i)*deltaX; 
x=(j-1)*deltaX; 
m=sym('m'); 
g=(((-1)^(m+1)+1)/m)*sin(m*pi*x/L1)*sinh(m*pi*y/L1)/sinh(m*pi*L2/L1); h=symsum(g,m,1,100); 
tx(i,j)=2*h*(tw2-tw1)/pi+tw1; 
end 
end 
tx 
迭代次数k = 695 
数值解温度场
数值解每次迭代的最大误差max1 =9.8243e-07 
解析温度场 tx = 
解析温度场
行的解析解和数值解的点
取第11行的解析解和数值解的点
行的解析解的直线，散点为其数值解的点 曲线为第11行的解析解的直线，散点为其数值解的点
解析解 第11行的误差=[数值解(11行) –解析解(11行)]/解析解
图像： 数值温度场 图像
图像 ：解析温度场tx图像：
数值解与解析解的误差
数值解与解析解的误差
程序内容： 数值计算matlab程序内容：>> t0=90; 
=10; 
L=10; 
c=0.25; 
p2=20; 
p1=p2/c; 
B=c*L; 
d=0.5*B; 
h=10; 
a=p2+1; 
b=p1+1; 
deltaX=B/p2; 
lambda=160; 
Bi=h*deltaX/lambda; 
=ones(a,b)*10; 
m1=ones(a,b)*3; 
m1(2:a-1,1)=zeros(a-2,1); 
m1(a,2:b-1)=ones(1,b-2); 
m1(1,2:b-1)=ones(1,b-2)*6; 
m1(2:a-1,b)=ones(a-2,1)*2; 
m1(1,b)=ones(1,1)*4; 
m1(a,b)=ones(1,1)*5; 
m1(1,1)=7; 
m1(a,1)=8; 
tn= ; 
max1=1.0; 
k=0; 
while ( max1>1e-6) 
k=k+1; 
max1=0; 
for i=1:1:a 
for j=1:1:b 
m=m1(i,j); 
n=tn(i,j); 
switch m 
case 0 
tn(i,j)=t0; 
case 1 
tn(i,j)=(2*tn(i-1,j)+tn(i,j-1)+tn(i,j+1)-4* )/(4+2*Bi)+ ; 
case 2 
tn(i,j)=(2*tn(i,j-1)+tn(i-1,j)+tn(i+1,j)-4* )/(4+2*Bi)+ ; 
case 3 
tn(i,j)=0.25*(tn(i,j-1)+tn(i,j+1)+tn(i-1,j)+tn(i+1,j)); 
case 4 
tn(i,j)=(tn(i,j-1)+tn(i+1,j)-2* )/(2*Bi+2)+ ; 
case 5 
tn(i,j)=(tn(i,j-1)+tn(i-1,j)-2* )/(2*Bi+2)+ ; 
case 6 
tn(i,j)=(2*tn(i+1,j)+tn(i,j-1)+tn(i,j+1)-4* )/(4+2*Bi)+ ; 
case 7 
tn(i,j)=t0; 
case 8 
tn(i,j)=t0; 
end 
er=abs(tn(i,j)-n); 
if er>max1 
max1=er; 
end 
end 
end 
=tn; 
end 
k 
ta=ones(a,b); 
Bi1=h*d/lambda; 
sbi=sqrt(Bi1); 
for i=1:1:a 
for j=1:1:b 
if i>(a+1)/2 
y=-(i-(a+1)/2)*deltaX; 
else y=((a+1)/2-i)*deltaX; 
end 
x=deltaX*(j-1); 
ta(i,j)=(cosh(sbi*(L-x)/d)+sbi*sinh(sbi*(L-x)/d))*(t0- )/(cosh(sbi*L/d)+sbi*sinh(sbi*L/d))+ ; 
end 
end 
ta 
迭代次数k =1461 
数值解温度场
解析温度场 ta 
解析温度场
行的解析解和数值解的点
取第11行的解析解和数值解的点
曲线为第11行的解析解的直线，散点为其数值解的点
行的解析解的直线，散点为其数值解的点
解析解 第11行的误差=[数值解(11行) –解析解(11行)]/解析解
图像如下
图像如下
数值温度场图像数值温度场图像 解析温度场图像解析温度场图像
数值解与解析解的误差数值解与解析解的误差
程序内容：
数值计算matlab程序内容：
>> tw=10; 
L2=15; 
c=0.75; 
L1=L2/c; 
p2=24 ; 
p1=p2/c; 
deltaX=2*L2/p2; 
a=p2+1; 
b=p1+1; 
lambda=16; 
qv0=24; 
=ones(a,b)*5; 
m1=ones(a,b); 
m1(1,:)=zeros(1,b); 
m1(2:a,b)=zeros(a-1,1); 
m1(2:a,1)=zeros(a-1,1); 
m1(a,2:b-1)=zeros(1,b-2); 
tn= ; 
max1=1.0; 
k=0; 
while(max1>1e-6) 
max1=0; 
k=k+1; 
for i=1:1:a 
for j=1:1:b 
m=m1(i,j); 
n=tn(i,j); 
switch m 
case 0 
tn(i,j)=tw; 
case 1 
tn(i,j)=0.25*(tn(i-1,j)+tn(i+1,j)+tn(i,j-1)+tn(i,j+1)+qv0*(deltaX^2)/lambda); 
end 
er=abs(tn(i,j)-n); 
if er>max1 
max1=er; 
end 
end 
end 
=tn; 
end 
k; 
tx=ones(a,b); 
for i=1:1:a 
for j=1:1:b 
if i>(a+1)/2 
y=-(i-(a+1)/2)*deltaX; 
else 
y=((a+1)/2-i)*deltaX; 
end 
if j>(b+1)/2 
x=(j-(b+1)/2)*deltaX; 
else 
x=-((b+1)/2-j)*deltaX; 
end 
m=sym('m'); 
xi=(2*m-1)*pi/2; 
g=((-1)^m)/(xi^3)*(cosh(xi*y/L1)/cosh(xi*L2/L1))*cos(xi*x/L1); h=symsum(g,m,1,100); 
tx(i,j)=2*qv0*L1^2/lambda*h+qv0*(L1^2-x^2)/(2*lambda)+tw; end 
end 
tx 
数值温度场
解析温度场tx 
取第13行的解析解和数值解的点
行的解析解和数值解的点
行的解析解的直线，散点为其数值解的点
曲线为第13行的解析解的直线，散点为其数值解的点
解析解 第13行的误差=[数值解(13行) –解析解(13行)]/解析解
数值温度场图像数值温度场图像 解析温度场图像解析温度场图像
数值解与解析解的误差数值解与解析解的误差。