實驗6平板穩(wěn)態(tài)溫度的計算

1、實驗6 平板穩(wěn)態(tài)溫度的計算6.1實驗指導(dǎo)本實驗在簡單的平板熱傳導(dǎo)模型下，利用線性方程組計算平板穩(wěn)態(tài)溫度的分布。當(dāng)選取節(jié)點較多時，MATLAB軟件的應(yīng)用就更顯重要。最后利用MATLAB的繪圖命令繪制平板的溫度分布，形象地表示出計算結(jié)果。表6.1給出了與本實驗相關(guān)的MATLAB命令或函數(shù)。表6.1 與本實驗相關(guān)的MATLAB命令或函數(shù)命令功能說明位置mod(i,N)計算i除以N的余例6.2if mod(i,N)=0判斷i是否為N的整數(shù)倍例6.2mesh(T)T為一個矩陣，以T的行標(biāo)為x軸值，以T的列標(biāo)為y軸值，以T的元素值為z軸值，繪制三維圖形。例6.26.2實驗內(nèi)容例6.1 在鋼板熱傳導(dǎo)的研究中

2、，常常用節(jié)點溫度來描述鋼板溫度的分布。假設(shè)圖6.1中鋼板已經(jīng)達(dá)到穩(wěn)態(tài)溫度分布，上下、左右四個邊界的溫度值如圖所示，而表示鋼板內(nèi)部四個節(jié)點的溫度。若忽略垂直于該截面方向的熱交換，那么內(nèi)部某節(jié)點的溫度值可以近似地等于與它相鄰四個節(jié)點溫度的算術(shù)平均值，如。請計算該鋼板的溫度分布。圖6.1 鋼板的節(jié)點分布（4個節(jié)點）解：根據(jù)已知條件可以得到以下線性方程組：化簡為標(biāo)準(zhǔn)的矩陣形式如下：在MATLAB命令窗口輸入： A=4,-1,-1,0; -1,4,0,-1; -1,0,4,-1; 0,-1,-1,4;b=30; 50; 60; 80;U=rref(A,b)結(jié)果為：U = 1.0000 0 0 0 21.

3、2500 0 1.0000 0 0 26.2500 0 0 1.0000 0 28.7500 0 0 0 1.0000 33.7500得到方程組的解為：，。例6.2 在例6.1中，把鋼板內(nèi)部分成了2×2個節(jié)點，本例把鋼板內(nèi)部分為5×5個節(jié)點，如圖6.2所示。求鋼板的穩(wěn)態(tài)溫度分布，并繪制溫度分布圖形。圖6.2 鋼板的節(jié)點分布（25個節(jié)點）解：根據(jù)例6.1中的討論，知：5×5個節(jié)點就構(gòu)成了一個具有25個方程的線性方程組。在MATLAB的M文件編輯器中編寫la19.m文件： % 計算鋼板的穩(wěn)態(tài)溫度分布clear close allN=input('N='

4、); % 輸入節(jié)點數(shù)，共有N×N個節(jié)點t_u=input('temperature up:'); % 輸入四個邊界的溫度值t_d=input('temperature down:');t_l=input('temperature left:');t_r=input('temperature right:');A=zeros(N*N);b=zeros(N*N,1); % 構(gòu)造N2×N2零矩陣A；構(gòu)造N2維零向量bfor i=1:N*N % 矩陣A的主對角線元素都是4 A(i,i)=4;endfor i=1:N*N

5、; % 給矩陣A和向量b賦值 if i<=N % 給向量b中和上邊界節(jié)點對應(yīng)的分量賦值 b(i)=t_u; end if mod(i,N)=0 % 給向量b中和右邊界節(jié)點對應(yīng)的分量賦值 b(i)=b(i)+t_r; end if mod(i,N)=1 % 給向量b中和左邊界節(jié)點對應(yīng)的分量賦值 b(i)=b(i)+t_l; end if i>N*(N-1) % 給向量b中和下邊界節(jié)點對應(yīng)的分量賦值 b(i)=b(i)+t_d; end if i>N % 給矩陣A中和上邊界無關(guān)的節(jié)點所對應(yīng)的元素賦值 A(i,i-N)=-1; end if mod(i,N)=1 % 給矩陣A中和左

6、邊界無關(guān)的節(jié)點所對應(yīng)的元素賦值 A(i,i-1)=-1; end if mod(i,N)=0 % 給矩陣A中和右邊界無關(guān)的節(jié)點所對應(yīng)的元素賦值 A(i,i+1)=-1; end if i<=N*(N-1) % 給矩陣A中和下邊界無關(guān)的節(jié)點所對應(yīng)的元素賦值 A(i,i+N)=-1; endendU=rref(A,b); % 對增廣矩陣進(jìn)行行初等變換，化為最簡行階梯矩陣for i=1:N % 把矩陣U的最后一列按矩陣形式賦給Bfor j=1:NB(i,j)=U(N*(i-1)+j,N*N+1);endendT(2:N+1,2:N+1)=B; % 把鋼板內(nèi)部溫度和和四周溫度值放入矩陣T中T(1

7、,2:N+1)=t_u;T(N+2,2:N+1)=t_d;T(2:N+1,1)=t_l;T(2:N+1,N+2)=t_r;T(1,N+2,1,N+2)=NaN; % 矩陣四個角沒有溫度值，故把非數(shù)NaN放入T % 顯示計算結(jié)果mesh(T) % 對鋼板溫度值繪制成曲面圖形在MATLAB的命令窗口中輸入：la19人機對話及運算結(jié)果為：N=5temperature up:20temperature down:50temperature left:10temperature right:30T = NaN 20.0000 20.0000 20.0000 20.0000 20.0000 NaN 10.

8、0000 16.5652 19.8667 21.8400 23.3415 25.3125 30.0000 10.0000 16.3958 21.0606 24.1538 26.2121 27.9111 30.0000 10.0000 17.9583 23.8261 27.5000 29.4426 30.1194 30.0000 10.0000 21.6087 28.7879 32.5769 33.9394 33.1233 30.0000 10.0000 29.6875 37.1395 40.0833 40.6140 38.4348 30.0000 NaN 50.0000 50.0000 50.

9、0000 50.0000 50.0000 NaN鋼板的溫度分布如圖6.3所示。其中x、y坐標(biāo)分別表示鋼板橫、縱方向的節(jié)點數(shù)，高度表示節(jié)點的溫度值，該三維圖形形象地反映了鋼板的溫度分布。圖6.3 鋼板的溫度分布6.3 實驗習(xí)題1、如圖6.4所示，假設(shè)鋼板已經(jīng)達(dá)到穩(wěn)態(tài)溫度分布。分別用來表示鋼板的上、下、左、右四個邊界的溫度值。請編寫一個通用Matlab程序，求出鋼板內(nèi)部六個節(jié)點的溫度值，要求四個邊界的溫度值由用戶輸入。（1），單位；（2），單位；（3）單位。圖6.4 鋼板的節(jié)點分布（6個節(jié)點）解：先列寫對此6點的方程：寫成矩陣形式程序為A=4,-1,0,-1,0,0;-1,4,-1,0,-1,0;0,-1,4,0,0,-1;-1,0,0,4,-1,0;0,-1,0,-1,4,-1;0,0,-1,0,-1,4, T_u=input('T_u= ') ,T_d=input('T_d= '),T_l=input('T_u= ') ,T_r=input('T_r= '),B=T_l+T_u;T_u; T_r+T_u; T_