传热学数值计算Word文件下载.docx
《传热学数值计算Word文件下载.docx》由会员分享，可在线阅读，更多相关《传热学数值计算Word文件下载.docx（47页珍藏版）》请在冰豆网上搜索。
y（j+1）=y（j）+1
p+1
q+1
n（j,i）=tw1+tw2*sinh（pi*y（j）/l1）*sin（pi*x（i）/l1）/sinh（pi*l2/l1）
各网格点用解析式得到的温度值：
解析二维温度分布图像：
误差分析：
取x=21，即位于板长一半处，温度随y（宽度）的变化曲线。
c1（:
1）取自于数值解，c1（:
2）取自于解析解
1）c1（:
2）
40.000040.0000
43.310643.4164
46.646546.8538
50.031350.3335
53.488953.8771
57.043057.5062
60.717861.2434
64.537665.1117
68.527369.1350
72.712273.3381
77.118777.7470
81.773682.3888
86.705087.2922
91.942392.4875
97.516298.0068
103.4592103.8840
109.8058110.1555
116.5925116.8600
123.8586124.0388
131.6461131.7363
140.0000140.0000
误差曲线：
由相对误差公式：
d1=（c1（:
2）-c1（:
1））./c1（:
2）可得：
d1=
0
0.0024
0.0044
0.0060
0.0072
0.0081
0.0086
0.0088
0.0085
0.0075
0.0067
0.0059
0.0050
0.0041
0.0032
0.0023
0.0015
0.0007
结论：
数值解与解析解吻合很好。
就x=21这一列，相对误差较小且在1%以内，但数值解较解析解偏小，且在平板中心附近的网格点的数值解较平板边缘数值解的相对误差大。
tw1=50%三边温度
tw2=100%一边温度
l1=20%板长20厘米
l2=10%板宽10厘米
m=41%划分成40*20的网格
k=1
c=zeros（n,m）
c（n,1）=（tw1+tw2）/2
c（n,m）=（tw1+tw2）/2
c（1,1）=tw1
c（1,m）=tw1
（m-1）
c（1,i）=tw1
c（n,i）=tw2
（n-1）
while（abs（k-c（j,i））>
tw1=50
m=40
n=20
l1=20
l2=10
tx=ones（20,40）
m+1
n+1
x=（i-1）*0.5
y=（j-1）*0.5
k=sym（'
k'
）
d=（（（-1）^（k+1）+1）/k）*sin（k*pi*x/l1）*sinh（k*pi*y/l1）/sinh（k*pi*l2/l1）
h=symsum（d,k,1,200）
tx（j,i）=2*h*（tw2-tw1）/pi+tw1
1）取自于数值解c1（:
1）=c1（:
2）=
50.000050.0000
51.979952.0881
53.973154.1851
55.990656.2998
58.043558.4409
60.141860.6165
62.295162.8344
64.511665.1016
66.798767.4243
69.162269.8077
71.606472.2558
74.133774.7710
76.744777.3546
79.437780.0056
82.209082.7214
85.052685.4977
87.960288.3278
90.921491.2035
93.924494.1151
96.955497.0513
100.000099.8408
误差分析曲线：
d2=abs（c1（:
d2=
0.0021
0.0039
0.0055
0.0068
0.0078
0.0091
0.0093
0.0092
0.0090
0.0079
0.0071
0.0062
0.0052
0.0042
0.0031
0.0020
0.0010
0.0016
误差结论：
就x=21这一列，相对误差较小且在1%以内，但数值解较解析解普遍偏小（最后一个数除外），且在平板中心附近的网格点的数值解较平板边缘数值解的相对误差大。
t0=200;
%肋基温度为200度
tf=0;
%环境温度为0度
L1=20;
%肋片长20厘米
L2=2;
%肋片高2厘米
h=0.01;
%对流换热系数单位：
w/（cm2*K）
a=11;
b=101;
%11*101的网格
dx=L1/（b-1）;
k=2;
%导热系数单位：
w/（cm*K）
Bi=h*dx/k;
ti=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;
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=ti;
max1=1.0;
w=0;
while（max1>
1e-6）
w=w+1;
max1=0;
fori=1:
a
forj=1:
b
m=m1（i,j）;
n=tn（i,j）;
switchm
case0
tn（i,j）=t0;
case1
tn（i,j）=（2*tn（i-1,j）+tn（i,j-1）+tn（i,j+1）-4*tf）/（4+2*Bi）+tf;
case2
tn（i,j）=（2*tn（i,j-1）+tn（i-1,j）+tn（i+1,j）-4*tf）/（4+2*Bi）+tf;
case3
tn（i,j）=0.25*（tn（i,j-1）+tn（i,j+1）+tn（i-1,j）+tn（i+1,j））;
case4
tn（i,j）=（tn（i,j-1）+tn（i+1,j）-2*tf）/（2*Bi+2）+tf;
case5
tn（i,j）=（tn（i,j-1）+tn（i-1,j）-2*tf）/（2*Bi+2）+tf;
case6
tn（i,j）=（2*tn（i+1,j）+tn（i,j-1）+tn（i,j+1）-4*tf）/（4+2*Bi）+tf;
case7
case8
end
er=abs（tn（i,j）-n）;
ifer>
max1
max1=er;
ti=tn;
ti
%假设垂直于纸面方向肋宽为1cm
q1=0;
qi=（tn（i,b）-tf）*h*0.2;
q1=q1+qi;
q1%x=20cm处肋片的散热量
q2=0;
b-1
q2j=（tn（a,j）-tf）*h*0.2;
q2=q2+q2j;
q2=q2*2;
q2%y=1cm与y=-1cm处肋片的散热量
q=q1+q2;
数值解肋片换热量：
q1=1.9022w/cm肋端沿y方向
q2=49.4938w/cm肋端沿x方向
q=51.3960w/cmq=q1+q2
Bi=h*1/k;
ta=ones（a,b）;
1:
ifi>
（a+1）/2
y=-（i-（a+1）/2）*dx;
elsey=（（a+1）/2-i）*dx;
x=dx*（j-1）;
ta（i,j）=（cosh（Bi^0.5*（L1-x）/1）+Bi^0.5*sinh（Bi^0.5*（L1-x）/1））*（t0-tf）/（cosh（Bi^0.5*L1/1）+Bi^0.5*sinh（Bi^0.5*L1/1））+tf;
ta
qi=（ta（i,b）-tf）*h*0.2;
q2j=（ta（a,j）-tf）*h*0.2;
解析解二维温度图像：
解析解肋片换热量：
q1=1.9006w/cm肋端沿y方向
q2=49.5481w/cm肋端沿x方向
q=51.4488w/cmq=q1+q2
取y=0，即位于板宽一半处，温度随x（长度）的变化曲线。
b1（:
1）取自于数值解b1（:
1）=b1（:
200.0000200.0000
197.5383197.4675
195.0883194.9745
192.6596192.5205
190.2587190.1051
187.8895187.7276
185.5543185.3877
183.2542183.0848
180.9896180.8186
178.7605178.5885
176.5667176.3942
174.4081174.2351
172.2842172.1109
170.1946170.0211
168.1389167.9654
166.1169165.9432
164.1280163.9542
162.1718161.9979
160.2480160.0741
158.3563158.1823
156.4961156.3222
154.6672154.4933
152.8692152.6953
151.1017150.9278
149.3644149.1906
147.6569147.4831
145.9788145.8052
144.3300144.1564
142.7099142.5365
141.1184140.9451
139.5550139.3818
138.0195137.8464
136.5115136.3386
135.0308134.8581
133.5771133.4046
132.1500131.9777
130.7494130.5772
129.3748129.2028
128.0261127.8543
126.7029126.5314
125.4051125.2337
124.1322123.9611
122.8842122.7133
121.6607121.4901
120.4615120.2911
119.2863119.1162
118.1350117.9651
117.0073116.8376
115.9029115.7335
114.8216114.6525
113.7633113.5945
112.7277112.5592
111.7147111.5464
110.7239110.5559
109.7552109.5875
108.8085108.6410
107.8835107.7163
106.9800106.8131
106.0979105.9312
105.2369105.0706
104.3970104.2309
103.5780103.4122
102.7796102.6140
102.0017101.8365
101.2442101.0793
100.5069100.3423
99.789799.6253
99.092598.9283
98.415098.2511
97.757197.5935
97.118896.9555
96.499996.3368
95.900295.7374
95.319795.1572
94.758394.5960
94.215794.0537
93.692093.5302
93.187093.0254
92.700692.5393
92.232792.0716
91.783291.6224
91.352091.1915
90.939190.7788
90.544390.3842
90.167790.0078
89.809089.6494
89.468389.3088
89.145488.9862
88.840388.6813
88.553088.3942
88.283488.1248
88.031487.8729
87.796987.6387
87.580087.4220
87.380687.2228
87.198787.0410
87.034186.8766
86.886986.7296
86.757186.6000
86.644686.4876
86.549486.3926
曲线适当放大：
d2=abs（b1（:
2）-b1（:
1））./b1（:
0.0004
0.0006
0.0008
0.0009
0.0011
0.0012
0.0013
0.0014
0.0016
0.0017
0.0018
0.0018
肋片总散热量数值解与解析解的相对误差为
q（解析解）=51.4488w/cm
q（数值解）=51.3960w/cm
t=abs（q（解析解）-q（数值解））/q（解析解）=0.001026
就y=0这行，相对误差较小且在0.2%以内，但数值解较解析解普遍偏大，且离肋基较近的网格点的数值解比离肋基较远的数值解相对误差大。
解析解的肋片总散热量与数值解吻合较好，但总散热量比数值解大。
tw=0;
%环境温度0度
L2=10;
%板宽的一半10厘米
c=0.5;
L1=L2/c;
%板长的一半20厘米
p2=40;
%划分成81*41的网格
p1=p2/c;
dx=2*L2/p2;
a=p2+1;
b=p1+1;
%导热系数单位w/（m2*K）
q=3;
%单位体积的产热率w/m3
ti=ones（a,b）*5;
m1=ones（a,b）;
m1（1,:
）=zeros（1,b）;
a,b）=zeros（a-1,1）;
a,1）=zeros（a-1,1）;
b-1）=zeros（1,b-2）;
y=0;
while（max1>
max1=0;
y=y+1;
m=m1（i,j）;
n=tn（i,j）;
switchm
case0
tn（i,j）=tw;
case1
tn（i,j）=0.25*（tn（i-1,j）+tn（i+1,j）+tn（i,j-1）+tn（i,j+1）+q*（dx^2）/k）;
er=abs（tn（i,j）-n）;
ifer>
ti=tn;
a1=ti（21:
41,41:
81）
tw1=0;
m=41;
n=21;
dx=L1/（m-1）;
tx=ones（n,m）;
m
x=（i-1）*0.5;
y=（j-1）*0.5;
p=sym（'
p'
）;
xi=（2*p-1）*pi/2;
d=（（-1）^p）/（xi^3）*（cosh（xi*y/L1）/cosh（xi*L2/L1））*cos（xi*x/L1）;
h=symsum（d,p,1,100）;
tx（j,i）=2*q*L1^2/k*h+q*（L1^2-x^2）/（2*k）+tw1;
tx
取x=1，即位于板长一半处，温度随y（宽度）的变化曲线。
68.312168.3231
68.145268.1562
67.644467.6553
66.809366.8200
65.6393
升级会员 