陳漓 莫小梅 孫美娟

(百色學院 a.材料科學與工程學院； b.化學與環(huán)境工程學院， 廣西 百色 533000)

?

運用Excel計算純流體飽和熱力學性質

陳漓a莫小梅b孫美娟a

(百色學院 a.材料科學與工程學院； b.化學與環(huán)境工程學院， 廣西 百色 533000)

利用Excel的數(shù)組公式和相關數(shù)據處理功能，通過簡單的設計，解決立方型狀態(tài)方程計算流體飽和熱力學性質所遇到的非線性高次方程、多次迭代問題，有效避開手工計算難度大和計算機計算需編程問題等，而且計算過程直觀、可視，利于掌握，計算結果準確、客觀。該方法可為研究者、學習者和教學者提供參考。

Excel；飽和蒸汽壓；改進PR狀態(tài)方程；熱力學性質

0 引言

計算純流體飽和熱力學性質，常常需要狀態(tài)方程，由于立方型狀態(tài)方程參數(shù)少且有較好的計算精度，因此得到廣泛的運用。在眾多立方形方程中PR狀態(tài)方程應用最為普遍。PR方程在預測液體體積上較之之前方程有所改善，但要進一步外推就比較困難。DAKEWU和SHULINCHEN修正PR方程，使方程的計算精度有顯著改進，在計算純流體的飽和蒸氣壓、飽和液相密度方面有良好的準確度。在此我們以改進PR狀態(tài)方程為模型結合Excel電子表格計算純流體的飽和熱力學性質。

純流體處于飽和狀態(tài)時，一般用4個狀態(tài)參量描述，溫度T、壓強p、汽相摩爾體積Vν和液相摩爾體積Vl。然而系統(tǒng)平衡時只有一個獨立變量，通常取溫度T或壓強p，故有兩類計算過程:第一類是蒸汽壓計算；第二類是沸點計算。

以第一類的蒸氣壓計算為例說明，并以改進PR狀態(tài)方程為模型[1]

(1)

純流體在飽和狀態(tài)時需滿足汽液平衡準則式:φν=φl

(2)

其中φν為汽相逸度，φl為液相逸度。

結合狀態(tài)方程，其逸度系數(shù)表達式可寫為：

(3)

純流體的汽液平衡的計算，首先要根據方程輸入臨界參數(shù)和偏心因子，計算給定溫度T下的改進PR方程常數(shù)a，b。假設p的初值，求狀態(tài)方程得到汽相摩爾體積Vν和液相摩爾體積Vl，由此判別方程式(2)是否滿足收斂條件，若不滿足，通過調節(jié)p，直到方程式(2)收斂，此時的p、Vν和Vl就是方程組式(1)和式(2)的解。由于在計算中涉及到高次方程運算，需借助計算機進行編程計算，操作過程較為繁瑣。如果能用辦公軟件Excel來解決，將給計算帶來很大的便利。

Excel是微軟辦公軟件的一個重要部分，它可以進行各種數(shù)據處理、統(tǒng)計分析等操作。對于熱力學數(shù)據的處理和運算也具有良好的交互界面，可簡單明了地展示計算原理與過程，有利于解決繁瑣的熱力學系統(tǒng)計算問題[2-5]。

1 純流體汽液平衡的理論基礎

1.1 計算飽和蒸汽壓

飽和蒸汽壓p進行首輪迭代計算時，要預先估計蒸汽壓的初值。運用Ednister偏心因子估算式ω=[3Tbr/7(1-Tbr)]lgpc-1能得到臨界參數(shù)和偏心因子估計蒸汽壓的初值的計算式[6]：

(4)

對于飽和蒸汽壓p的迭代式，可從式(2)利用Newton迭代法得到

p(k+1)=p(k)-[RT(lnφν-lnφ1)/(Vν-Vl)](k)

(5)

1.2 計算汽相、液相摩爾體積

對于立方形狀態(tài)方程摩爾體積的計算，當T

Vk+1=Vk-f(Vk)/f′(Vk)

(6)

為了方便迭代，將改進PR方程用摩爾體積的三次展開式來表示，即

(7)

f ′(V)為f(V)的一階導。

汽相摩爾體積Vν的求取通常以理想氣體體積V=RT/p為初值，而液相摩爾體積以V=b為初值，代入式(6)，得到Vl值后再代入等式右邊，一直迭代下去，直到滿足Vl<ε。利用Excel的計算功能，可方便同時迭代出Vν和Vl，這樣就極大提高計算的效率。

1.3 計算偏離焓、偏離熵等熱力學性質

流體的各種熱力學性質，比如焓、熵值和熱容等，可以利用偏離函數(shù)并結合狀態(tài)方程計算出來[7]。

(8)

(9)

(10)

偏離定壓摩爾熱容：

(11)

2 在Excel中進行純流體飽和熱力學性質的計算

為計算純流體汽液平衡下的飽和熱力學性質，在Excel中進行如下操作：

(1)以乙烷為例，分別在B3到E3單元格輸入臨界溫度、臨界壓強、氣體常數(shù)和偏心因子等數(shù)值。A6到F6單元格為方程常數(shù)a，b的值以及相應一階導數(shù)、二階導數(shù)的數(shù)值。

(2)在B10單元格輸入某一溫度“220”如圖1所示，即假設此時沸點溫度為220 K，運用式(4)計算飽和蒸氣壓的初值，即在A12單元格輸入“=C3*10^(7*(1+E3)/3*(1-B3/B10))”，把計算結果以粘貼數(shù)值的方式粘貼到A10單元格里。

(3)在單元格C9和D9分別輸入Vν的初值(V=RT/p)“=D3*B10/(A10*10^6)”和Vl的初值(V=b)“=B6”，單元格C10-C16和D10-D16分別輸入Newton迭代式(6)，如C10-C16單元格輸入如下數(shù)組公式“{=C9:C15-(C9:C15^3-(D3*B10/(A10*10^6)-0.645*B6)*C9:C15^2+1/(A10*10^6)*(D6-1.645*B6*D3*B10-2.29*A10*10^6*B6^2)*C9:C15-D6*B6/(A10*10^6)+0.645*D3*B10*B6^2/(A10*10^6)+0.645*B6^3)/(3*C9:C15^2-2*(D3*B10/(A10*10^6)-0.645*B6)*C9:C15+1/(A10*10^6)*(D6-1.645*B6*D3*B10-2.29*A10*10^6*B6^2))}”。同樣D10-D16單元格中式子與C10-C16內容相類似，所不同之處僅僅是把式中“C9:C15”改為“D9:D15”即可。一般來說大多數(shù)情況下迭代5-6次即可得到滿意結果,如圖1所示。為確保迭代的精度我們進行了7次迭代，單元格C16和D16得到的數(shù)值分別為汽相摩爾體積Vν和液相摩爾體積Vl，參見圖1。單元格F9和G9分別為式(7)f(V)的Vν和Vl收斂情況，f(V)應趨于零。

(4)接下來進一步計算逸度系數(shù)。在單元格F11和G11分別輸入計算汽相和液相逸度的關系式(3)，在單元格F11輸入“=A10*10^6*F10/(D3*B10)-1-LN(A10*10^6*(F10-B6)/(D3*B10))-D6/(2.299135707*B6*D3*B10)*LN((F10+1.9720678535*B6)/(F10-0.3270678535*B6))”。同理在單元格G11輸入液相逸度的關系式，計算后如果不滿足|lnφν-lnφl|≤ε則需要調整蒸汽壓p進行迭代計算，迭代關系式參見式(5)，即在單元格A14輸入式(5)，即“=(A10*10^6*(1-(F11-G11)/(A10*10^6*F10/(D3*B10)-A10*10^6*G10/(D3*B10))))/10^6”。把A14計算結果以粘貼數(shù)值的方式到A10單元格里，重復粘貼數(shù)值，Excel將對整個過程進行自動迭代計算，直至滿足|lnφν-lnφl|≤ε，一般迭代3-4次即能得到滿意的結果，此時得到p、Vν和Vl分別為飽和蒸汽壓、汽相摩爾體積和液相摩爾體積的數(shù)值。

(5)計算其它的熱力學性質。在單元格F13到F16和G13到G16分別輸入式(8)、式(9)、式(10)和式(11)，這樣我們同時得到乙烷汽相和液相偏離焓、偏離熵、偏離定容摩爾熱容和偏離定壓摩爾熱容等數(shù)值。例如，在單元格F13輸入“=A10*10^6*F10-D3*B10-1/(2.299135707*B6)*(D6-B10*E6)*LN((F10+1.9720678535*B6)/(F10-0.3270678535*B6))”，計算得到的數(shù)值為汽相偏離焓。

3 應用分析

為了檢驗Excel在計算純流體飽和性質的可靠性，我們以乙烷為研究對象，對該氣體飽和蒸汽壓等性質進行計算。在單元格B4、B5分別輸入乙烷的臨界溫度Tc和臨界壓強pc，單元格B7輸入偏心因子ω的數(shù)值參見圖1，計算乙烷溫度為160～305 K的飽和蒸汽壓、摩爾體積和偏離性質。計算的結果參見表1。

圖1 乙烷飽和性質的計算

溫度T/K飽和汽壓p/MPa摩爾體積汽相Vν/m3·mol-1液相Vl/m3·mol-1偏離焓汽相[HR]ν/J·mol-1液相[HR]l/J·mol-1偏離熵汽相[SR]ν/J·(mol·K)-1液相[SR]l/J·(mol·K)-1160002126225×10-25022×10-5-289508-156550863-01076-977709180007781875×10-25238×10-5-925438-150629562-03116-834806200021567300×10-35511×10-5-2285760-144285276-07085-717083220049043363×10-35867×10-5-4753996-137265643-13765-616090240096791728×10-36356×10-5-8817424-129190432-24176-525730260172119451×10-47081×10-5-15216885-119398562-40066-440765280282805238×10-48317×10-5-25472798-106428417-65514-354641300437522585×10-41153×10-4-45660183-84250988-118712-247348305484071903×10-41503×10-4-57951513-69336618-154363-191691

從表1計算的數(shù)據可以看出隨著溫度的增大，蒸發(fā)加劇密度增大飽和蒸汽壓升高，當溫度升高到305 K時，飽和汽壓逐步接近臨界壓強。同時溫度升高汽相摩爾體積Vν將隨著飽和汽壓升高而減小，液相摩爾體積Vl則隨著溫度升高而增大，在臨界點附近它們數(shù)值趨于一致。同樣偏離理想氣體性質焓和熵，汽相隨飽和汽壓升高，越遠離理想氣體狀態(tài)，液相隨溫度升高，偏離程度減小，在臨界點附近它們數(shù)值趨于一致。

把表1計算的飽和汽壓、摩爾體積、偏離焓和偏離熵結果與國內外計算熱力學性質軟件得到的數(shù)值是一致的，這也說明該方程能夠較好的描述純流體的熱力學性質。

上面的飽和熱力學性質的計算是以溫度T為獨立變量進行的，同樣也可以以蒸汽壓p為獨立變量來進行飽和熱力學性質的計算，只需在單元格B14輸入如下的迭代式即可

在單元格B14輸入“=B10*(1+(F11-G11)/((F13-G13)/(D3*B10)))”，多次迭代直到滿足|lnφν-lnφl|≤ε，此時溫度T為飽和蒸汽壓沸點的數(shù)值。

4 結論

通過上面的例子可以看出，運用Excel進行純流體汽液兩相飽和性質的計算，無需編程，只需在表格里輸入純流體的臨界性質和偏心因子的數(shù)值，通過簡單的迭代設計，Excel即可計算各種純流體的飽和熱力學性質。整個運算過程有良好的交互界面，并能有效地提高計算結果的準確度和可靠性。

Exce1具有強大的運算功能，但以往多用于簡單計算。實際上，充分利用Exce1提供的運算功能，只需輸入相應的關系式，就能完成復雜的熱力學工程計算，計算過程簡單、直觀，學生非常容易掌握。由于Exce1具有較強的數(shù)組公式和公式的復制功能，所以使用Exce1計算純流體熱力學性質具有無可比擬的優(yōu)越性。

[1]DAKEWU,SHULINCHEN.AmodifiedPeng-Robinsonequationonstate[J].ChemicalEngineeringCommunications， 1997, 156(1)：215-225.

[2]陳方.用EXCEL2000解決方程求根問題[J].湖南工業(yè)職業(yè)技術學院學報，2003，4(1)：14-16.

[3]班玉鳳，常圣泉，朱海峰，等.EXCEL在非理想系泡露點計算中的應用[J].計算機應用與軟件，2011，28(10)：275-277.

[4]于志家，陳傳棋，李香琴，等.應用Excel進行泡點與露點計算[J].化工高等教育，2012，126(4)：73-76.

[5]王智娟，胡粉娥，楊曉麗，等.Excel在化工熱力學計算中的應用[J].廣州化工，2013，41(3)：167-169.

[6]陳鐘秀，顧飛燕，胡望明，等.化工熱力學[M].3版.北京：化學工業(yè)出版社，2012：287-288.

[7]陳則韶.高等工程熱力學[M].2版.合肥：中國科學技術大學出版社，2014：120-125.

[責任編輯 劉景平]

Using Excel to Calculate the Thermodynamic Properties of Pure Fluid Saturation

CHEN Lia， MO Xiaomeib， SUN MeiJuana

(a.School of Materials Science and Engineering;b.School of Chemistry and Environmental Engineering, Baise University, Baise, Guangxi 533000, China)

The study uses Excel array formulas and data processing functions, through a simple design, to solve the problems of nonlinear equation of higher degree and multiple iterations faced when calculating thermodynamic properties of fluid saturation. It can effectively avoid the problems of manual calculation being difficult and computer calculation needing programming; the calculation process is intuitive, visual, conducive to grasp, and the calculation results are accurate and objective. The method provides a reference for researchers, learners and educators.

Excel; saturated vapor pressure; modified PR equation of state; thermodynamic properties

O642

A

1672-9021(2016)05-0113-05

陳漓(1962-)，男，廣西桂林人，百色學院材料科學與工程學院副教授，主要研究方向：熱力學。

廣西高?？茖W技術研究基金資助項目(2013YB244)。

2016-09-12