運(yùn)用Excel計(jì)算純流體飽和熱力學(xué)性質(zhì)

陳漓 莫小梅 孫美娟

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

?

運(yùn)用Excel計(jì)算純流體飽和熱力學(xué)性質(zhì)

陳漓a莫小梅b孫美娟a

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

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

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

0 引言

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

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

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

(1)

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

(2)

其中φν為汽相逸度，φl(shuí)為液相逸度。

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

(3)

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

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

1 純流體汽液平衡的理論基礎(chǔ)

1.1 計(jì)算飽和蒸汽壓

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

(4)

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

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

(5)

1.2 計(jì)算汽相、液相摩爾體積

對(duì)于立方形狀態(tài)方程摩爾體積的計(jì)算，當(dāng)T

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

(6)

為了方便迭代，將改進(jìn)PR方程用摩爾體積的三次展開式來(lái)表示，即

(7)

f ′(V)為f(V)的一階導(dǎo)。

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

1.3 計(jì)算偏離焓、偏離熵等熱力學(xué)性質(zhì)

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

(8)

(9)

(10)

偏離定壓摩爾熱容：

(11)

2 在Excel中進(jìn)行純流體飽和熱力學(xué)性質(zhì)的計(jì)算

為計(jì)算純流體汽液平衡下的飽和熱力學(xué)性質(zhì)，在Excel中進(jìn)行如下操作：

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

(2)在B10單元格輸入某一溫度“220”如圖1所示，即假設(shè)此時(shí)沸點(diǎn)溫度為220 K，運(yùn)用式(4)計(jì)算飽和蒸氣壓的初值，即在A12單元格輸入“=C3*10^(7*(1+E3)/3*(1-B3/B10))”，把計(jì)算結(jié)果以粘貼數(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內(nèi)容相類似，所不同之處僅僅是把式中“C9:C15”改為“D9:D15”即可。一般來(lái)說(shuō)大多數(shù)情況下迭代5-6次即可得到滿意結(jié)果,如圖1所示。為確保迭代的精度我們進(jìn)行了7次迭代，單元格C16和D16得到的數(shù)值分別為汽相摩爾體積Vν和液相摩爾體積Vl，參見圖1。單元格F9和G9分別為式(7)f(V)的Vν和Vl收斂情況，f(V)應(yīng)趨于零。

(4)接下來(lái)進(jìn)一步計(jì)算逸度系數(shù)。在單元格F11和G11分別輸入計(jì)算汽相和液相逸度的關(guān)系式(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輸入液相逸度的關(guān)系式，計(jì)算后如果不滿足|lnφν-lnφl(shuí)|≤ε則需要調(diào)整蒸汽壓p進(jìn)行迭代計(jì)算，迭代關(guān)系式參見式(5)，即在單元格A14輸入式(5)，即“=(A10*10^6*(1-(F11-G11)/(A10*10^6*F10/(D3*B10)-A10*10^6*G10/(D3*B10))))/10^6”。把A14計(jì)算結(jié)果以粘貼數(shù)值的方式到A10單元格里，重復(fù)粘貼數(shù)值，Excel將對(duì)整個(gè)過(guò)程進(jìn)行自動(dòng)迭代計(jì)算，直至滿足|lnφν-lnφl(shuí)|≤ε，一般迭代3-4次即能得到滿意的結(jié)果，此時(shí)得到p、Vν和Vl分別為飽和蒸汽壓、汽相摩爾體積和液相摩爾體積的數(shù)值。

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

3 應(yīng)用分析

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

圖1 乙烷飽和性質(zhì)的計(jì)算

溫度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計(jì)算的數(shù)據(jù)可以看出隨著溫度的增大，蒸發(fā)加劇密度增大飽和蒸汽壓升高，當(dāng)溫度升高到305 K時(shí)，飽和汽壓逐步接近臨界壓強(qiáng)。同時(shí)溫度升高汽相摩爾體積Vν將隨著飽和汽壓升高而減小，液相摩爾體積Vl則隨著溫度升高而增大，在臨界點(diǎn)附近它們數(shù)值趨于一致。同樣偏離理想氣體性質(zhì)焓和熵，汽相隨飽和汽壓升高，越遠(yuǎn)離理想氣體狀態(tài)，液相隨溫度升高，偏離程度減小，在臨界點(diǎn)附近它們數(shù)值趨于一致。

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

上面的飽和熱力學(xué)性質(zhì)的計(jì)算是以溫度T為獨(dú)立變量進(jìn)行的，同樣也可以以蒸汽壓p為獨(dú)立變量來(lái)進(jìn)行飽和熱力學(xué)性質(zhì)的計(jì)算，只需在單元格B14輸入如下的迭代式即可

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

4 結(jié)論

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

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

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

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

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

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

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

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

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

[責(zé)任編輯 劉景平]

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-)，男，廣西桂林人，百色學(xué)院材料科學(xué)與工程學(xué)院副教授，主要研究方向：熱力學(xué)。

廣西高?？茖W(xué)技術(shù)研究基金資助項(xiàng)目(2013YB244)。

2016-09-12