传热学大作业.docx
《传热学大作业.docx》由会员分享，可在线阅读，更多相关《传热学大作业.docx（14页珍藏版）》请在冰点文库上搜索。
传热学大作业
传热学第四章大作业
二维稳态导热问题
的数值解法
姓名：
班级：
学号：
第一题：
如图所示，一个无限长矩形柱体，其横截面的边长分别为
和
，常物性。
该问题可视为二维稳态导热问题，边界条件如图中所示，其中
=0.6m，
=0.4m，
=60℃，
=20℃，
。
1）编写程序求解二维导热方程。
2）绘制x=
/2和y=
/2处的温度场，并与解析解进行比较。
已知矩形内的温度场的解析解为：
解：
（1）建立控制方程及定解条件
控制方程：
定解条件：
（2）区域离散化（确立节点）
将矩形区域分为M*N个网格，其中x方向上的步长
;y方向的步长
。
设节点为（m,n）。
（3）建立节点离散方程
对节点（m,n）有：
内节点:
化简得
边界节点：
（4）编程求解，程序见附录
取M=N=50得到矩形区域各节点温度（见附件一），为方便在这里仅给出M=N=10时温度分布数据，如下表：
60
60
60
60
60
60
60
60
60
60
60
60.44865
60.82879
61.09192
61.2134
61.18838
61.02741
60.75314
60.39773
60
60
60.92774
61.70958
62.24681
62.4919
62.43718
62.10517
61.54236
60.81427
60
60
61.47168
62.69908
63.53154
63.90366
63.80923
63.2855
62.40502
61.26912
60
60
62.12418
63.86412
65.02086
65.52276
65.37208
64.62439
63.3812
61.78318
60
60
62.94745
65.28836
66.80122
67.43191
67.20027
66.18351
64.51509
62.37953
60
60
64.0403
67.08063
68.97364
69.72431
69.3774
68.03224
65.8566
63.08431
60
60
65.57759
69.38281
71.65378
72.50452
71.9982
70.25008
67.46328
63.92776
60
60
67.90258
72.36689
74.96513
75.88788
75.17091
72.92955
69.40274
64.94553
60
66.18034
71.75571
76.18034
79.02113
80
79.02113
76.18034
71.75571
66.18034
60
M=N=50时矩形区域节点温度分布图如图：
M=N=50时x=
/2处数值解与解析解温度场分布如下：
M=N=10时将数值解与理论解进行比较，如下表所示：
y
0
0.4/9
0.4/9*2
0.4/9*3
0.4/9*4
0.4/9*5
0.4/9*6
0.4/9*7
0.4/9*8
0.4
数值解
60
61.2134
62.4919
63.90366
65.52276
67.43191
69.72431
72.50452
75.88788
80
理论解
60
61.17447
62.41282
63.78244
65.35781
67.22464
69.48449
72.26029
75.70303
80
相对误差
0
0.064%
0.127%
0.190%
0.252%
0.308%
0.345%
0.338%
0.244%
0
M=N=50时y=
/2处数值解与解析解温度场分布如下图：
M=N=10时将数值解与理论解进行比较，如下表所示：
x
0
0.6/9
0.6/9*2
0.6/9*3
0.6/9*4
0.6/9*5
0.6/9*6
0.6/9*7
0.6/9*8
0.6
数值解
60
62.12418
63.86412
65.02086
65.52276
65.37208
64.62439
63.3812
61.78318
60
理论解
60
62.13724
64.0167
65.41169
66.15395
66.15395
65.41169
64.0167
62.13724
60
相对误差
0
0.021%
0.238%
0.598%
0.954%
1.182%
1.204%
0.993%
0.569%
0
第二题：
将第一题中
处的边界条件变为
，其他条件不变。
1）编写程序求解二维导热方程并计算从y=0处导入的热量
。
2）当
时，该二维导热问题可简化为一维导热问题。
在一维的近似下，试计算从y=0处导入的热量
，并比较不同L2/L1下
的比值。
由该问题的解析解可知：
L2/L1
0.007
0.01
0.05
0.08
0.1
0.9987
0.9912
0.956
0.93
0.912
解：
编程得到M=N=50温度分布数据见附件二，得到温度分布图如下：
这里也仅给出M=N=10时温度分布数据，见下表：
60
60
60
60
60
60
60
60
60
60
60
58.92484
58.03075
57.41692
57.11051
57.11051
57.41692
58.03075
58.92484
60
60
57.76922
55.93693
54.6972
54.08483
54.08483
54.6972
55.93693
57.76922
60
60
56.43648
53.57975
51.69865
50.787
50.787
51.69865
53.57975
56.43648
60
60
54.78962
50.78896
48.26924
47.08399
47.08399
48.26924
50.78896
54.78962
60
60
52.6051
47.33997
44.24672
42.8542
42.8542
44.24672
47.33997
52.6051
60
60
49.47402
42.9257
39.46834
38.00551
38.00551
39.46834
42.9257
49.47402
60
60
44.57508
37.13768
33.80348
32.50667
32.50667
33.80348
37.13768
44.57508
60
60
36.12615
29.52601
27.23313
26.43148
26.43148
27.23313
29.52601
36.12615
60
20
20
20
20
20
20
20
20
20
20
从y=0处导入的热量
可以近似看作从y=0向y=0.4/9处传递的热量，y=0.4/9处的温度分别为：
x
0
0.6/9
0.6/9*2
0.6/9*3
0.6/9*4
0.6/9*5
0.6/9*6
0.6/9*7
0.6/9*8
0.6
T
60
58.92484
58.03075
57.41692
57.11051
57.11051
57.41692
58.03075
58.92484
60
则
整理得
当
时，该二维导热问题可简化为一维导热问题。
当
时，令M=N=100，得
=1129571.429
看作一维导热时
不同
下
比值如下表：
L2/L1
0.007
0.01
0.05
0.08
0.1
1129571.429
789540
152812
92890
72908
1142857.143
800000
160000
100000
80000
0.9884
0.9869
0.9551
0.9289
0.9114
0.9987
0.9912
0.956
0.93
0.912
附录
Matlab程序：
第一题：
functionchuanredazuoye
clear
clc
L1=0.6;%矩形长度
L2=0.4;%矩形宽度
Tw1=60;
Tw2=20;
globalMN%设置网格数
M=input（'请输入将区间[0,L1]等分的个数M:
'）;
N=input（'请输入将区间[0,L2]等分的个数N:
'）;
XDIF=L1/M;%x方向上的步长
YDIF=L2/N;%y方向上的步长
axis（[0,M,0,N]）;
grid%设置网格
U=initial（M,N,Tw1,Tw2,L1,XDIF）;
AIM=YDIF/XDIF;
AIP=YDIF/XDIF;
AJM=XDIF/YDIF;
AJP=XDIF/YDIF;
AP=2*YDIF/XDIF+2*XDIF/YDIF;
CON=0;%离散方程系数
T=rechuandao（U,M,N,AIM,AIP,AJM,AJP,AP）;
A=flipud（T）;
mesh（A）;
title（'温度分布'）
xlabel（'x'）
ylabel（'y'）
zlabel（'T'）
functionU=rechuandao（U,M,N,AIM,AIP,AJM,AJP,AP）
while1
temp=U;
fori=2:
M-1
forj=2:
N-1
U（i,j）=（AIM*U（i,j-1）+AIP*U（i,j+1）+AJM*U（i-1,j）+AJP*U（i+1,j））/AP;%设置迭代条件
end
end
eps=abs（U-temp）;
ifmax（max（eps））<1e-8
break;%限制迭代次数
end
end
functionU=initial（M,N,Tw1,Tw2,L1,XDIF）
U=zeros（M,N）;%赋温度矩阵初值
U（:
1）=Tw1;
U（:
M）=Tw1;
U（N,:
）=Tw1;
fori=1:
M
U（1,i）=Tw1+Tw2*sin（pi*i*XDIF/L1）;
end%初始和边界条件的设定
第二题：
functiongaibianjietiaojian
clear
clc
L1=0.6;%矩形长度
L2=0.4;%矩形宽度
Tw1=60;
Tw2=20;
globalMN%设置网格数
M=input（'请输入将区间[0,L1]等分的个数M:
'）;
N=input（'请输入将区间[0,L2]等分的个数N:
'）;
XDIF=L1/M;%x方向上的步长
YDIF=L2/N;%y方向上的步长
axis（[0,M,0,N]）;
grid%设置网格
U=initial（M,N,Tw1,Tw2）;
AIM=YDIF/XDIF;
AIP=YDIF/XDIF;
AJM=XDIF/YDIF;
AJP=XDIF/YDIF;
AP=2*YDIF/XDIF+2*XDIF/YDIF;
CON=0;%离散方程系数
T=rechuandao（U,M,N,AIM,AIP,AJM,AJP,AP）;
A=flipud（T）;
mesh（A）;
title（'改变边界条件后温度分布'）
xlabel（'x'）
ylabel（'y'）
zlabel（'T'）
functionU=rechuandao（U,M,N,AIM,AIP,AJM,AJP,AP）
while1
temp=U;
fori=2:
M-1
forj=2:
N-1
U（i,j）=（AIM*U（i,j-1）+AIP*U（i,j+1）+AJM*U（i-1,j）+AJP*U（i+1,j））/AP;
%设置迭代条件
end
end
eps=abs（U-temp）;
ifmax（max（eps））<1e-8
break;%限制迭代次数
end
end
functionU=initial（M,N,Tw1,Tw2）
U=zeros（M,N）;%赋温度矩阵初值
U（:
1）=Tw1;
U（:
M）=Tw1;
U（N,:
）=Tw1;
U（1,:
）=Tw2;%初始和边界条件的设定
升级会员 