傳熱學(xué)上機(jī)作業(yè)-墻角溫度場(chǎng)分布的數(shù)值模擬

《傳熱學(xué)》上機(jī)實(shí)踐大作業(yè)二維導(dǎo)熱物體溫度場(chǎng)的數(shù)值模擬能動(dòng)A02趙凱2010031134一、物理問題有一個(gè)用磚砌成的長(zhǎng)方形截面的冷空氣通道，其截面尺寸如下圖所示，假設(shè)在垂直于紙面方向上冷空氣及磚墻的溫度變化很小，可以近似地予以忽略。在下列兩種情況下試計(jì)算:磚墻橫截面上的溫度分布；垂直于紙面方向的每米長(zhǎng)度上通過磚墻的導(dǎo)熱量。第一種情況：內(nèi)外壁分別均勻地維持在0℃及30℃；第二種情況：內(nèi)外表面均為第三類邊界條件，且已知：t=30℃,/i=10.35W/m2.^col 1t=10℃,/i=3.93W/m^Koo2 2石專墻導(dǎo)熱系猴=0.35/加二、數(shù)學(xué)描寫1、控制方程該問題為無內(nèi)熱源的二維穩(wěn)態(tài)導(dǎo)熱問題，因此控制方程為導(dǎo)熱微分方程:d21 d21 人——十——=0dxx2 dyy22、邊界條件該問題中，導(dǎo)熱物體在x方向上，y方向上都是對(duì)稱的，因此可以只取其中的四分之一部分作為研究對(duì)象，其他部分情況完全相同，如下圖所示:對(duì)于上圖所示各邊界:邊界1：由對(duì)稱性可知：其為絕熱邊界，即邊界1：由對(duì)稱性可知：其為絕熱邊界，即q=0。w邊界2：第一種情況：其為等溫邊界，滿足第一類邊界條件。即:t=0。C即:w第二種情況：其為對(duì)流邊界，滿足第三類邊界條件。即:邊界3：邊界3：第一種情況：其為等溫邊界，滿足第一類邊界條件。即:第二種情況：其為對(duì)流邊界，滿足第三類邊界條件。第二種情況：其為對(duì)流邊界，滿足第三類邊界條件。,at即:q=—X(一)=h(t—t)

w anwiwf即:三、方程離散如下圖所示，用一系列與坐標(biāo)軸平行的間隔10cm的網(wǎng)格線將求解區(qū)域劃分成子區(qū)域。111可將上圖所示各節(jié)點(diǎn)分成內(nèi)節(jié)點(diǎn)與邊界點(diǎn)兩類。分別利用熱平衡法列各個(gè)節(jié)點(diǎn)的代數(shù)方程。第一種情況（等溫邊界）：邊界點(diǎn)：,1邊界1（絕熱邊界）：/1— （2t）+t-+t1J，m—2~5乂，【：八、、乂，【人m,1 4m,2 m_1,1 m+1,11_ 、…，t= （21 +1 +1 ）,n=8~1116,n4 15,n16,n-1 16,n+1邊界2（內(nèi)等溫邊界）：t=0,m=6,n=1~7;m=7~16,n=7m,n邊界3（外等溫邊界）：t=30,m=1,n―1~12;m=2~16,n=12m,n內(nèi)節(jié)點(diǎn):， 1,， ， ， ，、.t=（t +1 +1 +1 ）;m,n 4m+1,n m-1,n m,n+1 m,n-1m=2~5,n=2~11;m=6~15,n=8~11第二種情況（對(duì)流邊界）：邊界點(diǎn)：t邊界1（絕熱邊界）：m,11(21t邊界1（絕熱邊界）：m,11(214m,2+1m—1,1+1 ),m=2~5m+1,1t16,n1(2t4 15,n+t16,n—1+1 ),n=8~1116,n+12t+t+t+2Bit邊界2（內(nèi)對(duì)流邊界）：t6n5,n 6,n+16,n—1A邊界2（內(nèi)對(duì)流邊界）：t6n2(Bi+2)A2t+t+t +2Bitt=m,8 m+1,7 m—1,7 A11,m=7~16m,7 2(Bi+2)A12t+t+t+2Bitt= 2,n 1,n+1 1,n—1 A22.n=1~11邊界3(外對(duì)流邊界)：1,n 2(Bi+2) ,A22t +t+t +2Bitt=m,11 m+1,12 m—1,12 A22,m=2~16m,12 2(Bi+2)A22(t+t)+t+t+2Bitt— 5,7 6,8 7,7 6,6A11內(nèi)角點(diǎn)： 6,7 2(Bi+3)A1t+t+2Bitt―2,12 1,11 A22外角點(diǎn)： 1,12 2(Bi+1)A內(nèi)節(jié)點(diǎn)：(t +t+t+t);m+1,n m—1,n m,n+1 m,n—1m—2~5,n—2~11;m—6~15,n—8~11(BiA1,t—t—10；BiA281—30)四、編程求解第一種情況(等溫邊界)：Fortran程序代碼如下所示：Programdenwengimplicitnoneinteger::t1=0integer::t2=30integerm,nreal::t(16,12),ta(16,12),et(16,12)real::epslona=1realfainei,fainei1,fainei2,fainei3,fainei4,fainei5,fainei6,fainei7realfaiwai,faiwai1,faiwai2,faiwai3,faiwai4,faiwai5,faiwai6,faiwai7realpianchadon=1,7t(6,n)=t1enddodom=7,16t(m,7)=t1enddodon=1,12t(1,n)=t2enddodom=2,16t(m,12)=t2enddodom=2,5don=1,11t(m,n)=10enddoenddodom=6,16don=8,11t(m,n)=10enddoenddoopen(01,file='dengwen.dat')dowhile(epslona>0.00000001)dom=2,5ta(m,1)=0.25*(2*t(m,2)+t(m-1,1)+t(m+1,1))enddodom=2,5don=2,11ta(m,n)=0.25*(t(m+1,n)+t(m-1,n)+t(m,n+1)+t(m,n-1))enddoenddodom=6,15don=8,11ta(m,n)=0.25*(t(m+1,n)+t(m-1,n)+t(m,n+1)+t(m,n-1))enddoenddodon=8,11ta(16,n)=0.25*(2*t(15,n)+t(16,n-1)+t(16,n+1))enddodon=1,7ta(6,n)=t1enddodom=7,16ta(m,7)=t1enddodon=1,12ta(1,n)=t2enddodom=2,16ta(m,12)=t2enddodom=1,16don=1,12et(m,n)=abs(ta(m,n)-t(m,n))enddoenddoepslona=maxval(et(1:16,1:12))dom=1,16don=1,12t(m,n)=ta(m,n)enddoenddoenddofainei1=0.5*lanbuda*t(5,1)fainei3=lanbuda*t(5,8)fainei5=0.5*lanbuda*t(16,8)fainei2=0don=2,7fainei6=lanbuda*t(5,n)fainei2=fainei2+fainei6enddofainei4=0dom=6,15fainei7=lanbuda*t(m,8)fainei4=fainei4+fainei7enddofainei=4*(fainei1+fainei2+fainei3+fainei4+fainei5)faiwai1=0.5*lanbuda*(30-t(2,1))faiwai3=lanbuda*(30-t(2,11))faiwai5=0.5*lanbuda*(30-t(16,11))faiwai2=0don=2,10faiwai6=lanbuda*(30-t(2,n))faiwai2=faiwai2+faiwai6enddofaiwai4=0dom=3,15faiwai7=lanbuda*(30-t(m,11))faiwai4=faiwai4+faiwai7enddofaiwai=4*(faiwai1+faiwai2+faiwai3+faiwai4+faiwai5)print*,' m nt'dom=1,16don=1,12print*,m,n,t(m,n)write(01,*)m,n,t(m,n)enddoenddopiancha=abs(fainei-faiwai)/((fainei+faiwai)/2)print*,'內(nèi)部熱流量=',faineiprint*,'外部熱流量=',faiwaiprint*,'熱平衡偏差=',pianchaendprogramdenweng運(yùn)行結(jié)果如圖所示:'C:\ProgramFiles\Mi:crosDftVisualStudio\Comm□ SDEV98\M￥PROJECTS\Debug\zui3....0,產(chǎn)mnt1130.09S091230.03@031330.096091430.09S091620.036031630.096091730.0BS001a20.BHGBB1930.0960911030.0BQ0011130-BBQB811230.096Q3Z1Z4.1D2662224.121402324.1S3G32424.308452524.529982G24.S9G442725.456812826.222052927.ll.7tl21028.Q690621129.03^522123fl-0BGBH|*C:\ProgramFiles\MicrasoftVisualStudio\Ccirrmon\MSDEV9B\MVPROJECT￡\Debug\3.eKeb40.00O000@(E+005B.00B0QBB(E+0Q60.00O00O8(E+007B.OOOBOOBE+OO8 6.0109469 12.8176710 18.0176311 24.0108712 30.OQ00SI1O.00OO0OBIE+0020.000O000E+003B.00fl0BOB(E+QQ40.00O00OS(E1035B.EHtJEBtJBtE+En60.000000S(E+0370.00O000BIE+008 6.0091819 12.01483IQ1S.014S111 24.00914雪窗a雪窗a部部平ss內(nèi)are270,6282239.GGS10.1213418tocontinue第二種情況(對(duì)流邊界):Fortran程序代碼如下所示：programduiliuimplicitnoneinteger::t1=10integer::t2=30integerm,nreal::t(16,12),ta(16,12),et(16,12)real::epslona=1realbi1,bi2realfainei,fainei1,fainei2,fainei3,fainei4,fainei5,fainei6,fainei7realfaiwai,faiwai1,faiwai2,faiwai3,faiwai4,faiwai5,faiwai6,faiwai7realpianchabi1=h1*detax/lanbudabi2=h2*detax/lanbudadom=1,16don=1,12t(m,n)=10enddoenddoopen(01,file='crs.dat')dowhile(epslona>0.000000001)dom=2,5ta(m,1)=0.25*(2*t(m,2)+t(m-1,1)+t(m+1,1))enddodon=8,11ta(16,n)=0.25*(2*t(15,n)+t(16,n-1)+t(16,n+1))enddota(6,n)=(2*t(5,n)+t(6,n+1)+t(6,n-1)+2*bi1*t1)/(2*bi1+4)enddodom=7,15ta(m,7)=(2*t(m,8)+t(m+1,7)+t(m-1,7)+2*bi1*t1)/(2*bi1+4)enddodon=2,11ta(1,n)=(2*t(2,n)+t(1,n+1)+t(1,n-1)+2*bi2*t2)/(2*bi2+4)enddodom=2,15ta(m,12)=(2*t(m,11)+t(m+1,12)+t(m-1,12)+2*bi2*t2)/(2*bi2+4)enddodom=2,5don=2,11ta(m,n)=0.25*(t(m+1,n)+t(m-1,n)+t(m,n+1)+t(m,n-1))enddoenddodom=6,15don=8,11ta(m,n)=0.25*(t(m+1,n)+t(m-1,n)+t(m,n+1)+t(m,n-1))enddoenddodon=2,6ta(6,7)=(2*t(5,7)+2*t(6,8)+t(7,7)+t(6,6)+2*bi1*t1)/(2*bi1+6)ta(1,12)=(t(2,12)+t(1,11)+2*bi2*t2)/(2*bi2+2)ta(6,1)=(t(5,1)+t(6,2)+bi1*t1)/(bi1+2)ta(16,7)=(t(16,8)+t(15,7)+bi1*t1)/(bi1+2)ta(16,12)=(t(16,11)+t(15,12)+bi2*t2)/(bi2+2)ta(1,1)=(t(2,1)+t(1,2)+bi2*t2)/(bi2+2)dom=1,16don=1,12et(m,n)=abs(ta(m,n)-t(m,n))enddoenddoepslona=maxval(et(1:16,1:12))dom=1,16don=1,12t(m,n)=ta(m,n)enddoenddoenddofainei1=0.05*h1*(t(6,1)-10)fainei3=0.1*h1*(t(6,7)-10)fainei5=0.05*h1*(t(16,7)-10)fainei2=0don=2,6fainei6=0.1*h1*(t(6,n)-10)fainei2=fainei2+fainei6enddofainei4=0dom=7,15fainei7=0.05*h1*(t(m,8)-10)fainei4=fainei4+fainei7enddofainei=4*(fainei1+fainei2+fainei3+fainei4+fainei5)faiwai6=0.1*h2*(30-t(1,n))faiwai2=faiwai2+faiwai6enddofaiwai4=0dom=2,15faiwai7=0.1*h2*(30-t(m,12))faiwai4=faiwai4+faiwai7enddofaiwai=4*(faiwai1+faiwai2+faiwai3+faiwai4+faiwai5)don=1,12dom=1,16print*,m,n,t(m,n)write(01,*)m,n,t(m,n)enddoenddopiancha=abs(fainei-faiwai)/((fainei+faiwai)/2)print*,'內(nèi)部熱流量=',faineiprint*,'外部熱流量=',faiwaiprint*,'熱平衡偏差=',pianchaclose(01)endprogramduiliufaiwai1=0.05*h2*(30-t(1,1))faiwai3=0.1*h2*(30-t(1,12))faiwai5=0.05*h2*(30-t(16,12))faiwai2=0don=2,11運(yùn)行結(jié)果如圖所示:'■'C:\Pnognam rcrasoftVr&jalStudio\Comman\NSDEV9S\MYPROJECTS\Debug\zui4....mnt1128.586092125.817023123.0O2124120.118625117.15t206114.12212?1Q.0Q0BQ00E+0QS1Q-BBBBQQ8E+BQ916.0@00603￡+06101S.0S06@00E+0S111G-BBBBGBBE+BB1216.SS0iBS0dE-|-0S131B.OEJEJ0O0OE+OEJ141151S.0S0IB603E+0S1610.OS00S09E+O01228.594052225.839933223.036424220.158075217.192046214.14-464?