趙 丹 丁國良 胡海濤

(上海交通大學制冷與低溫工程研究所 上海 200240)

質(zhì)量和能量嚴格守恒的蒸發(fā)器動態(tài)仿真模型

趙 丹 丁國良 胡海濤

(上海交通大學制冷與低溫工程研究所 上海 200240)

現(xiàn)有用于蒸發(fā)器動態(tài)仿真的移動邊界模型都是基于質(zhì)量和能量守恒方程展開獲得，而展開過程中存在誤差，從而導致了現(xiàn)有模型質(zhì)量和能量不嚴格守恒, 進而影響仿真精度。為了保證求解過程質(zhì)量和能量的嚴格守恒，本文選取蒸發(fā)器內(nèi)制冷劑總質(zhì)量和總能量作為控制方程的狀態(tài)參數(shù)，直接將質(zhì)量和能量守恒方程作為控制方程求解，從而避免了方程展開的誤差。為了在已知制冷劑總質(zhì)量和總能量的情況下求解蒸發(fā)器內(nèi)制冷劑分布情況和制冷劑狀態(tài)，本文開發(fā)了制冷劑不同分布情況的計算公式和計算方法。仿真案例表明，新方法在仿真48 h蒸發(fā)器性能中，質(zhì)量和能量嚴格守恒，計算穩(wěn)定，仿真結(jié)果和實驗結(jié)果吻合。

蒸發(fā)器；動態(tài)仿真；移動邊界模型

蒸氣壓縮式制冷系統(tǒng)廣泛應用于冰箱和空調(diào)系統(tǒng)。蒸發(fā)器直接參與系統(tǒng)換熱環(huán)節(jié)，是制冷系統(tǒng)中最重要的部件之一，所以蒸發(fā)器結(jié)構(gòu)和控制策略設(shè)計對系統(tǒng)性能的影響至關(guān)重要。動態(tài)仿真方法是實現(xiàn)蒸發(fā)器結(jié)構(gòu)和控制策略設(shè)計的有效方法，相對于傳統(tǒng)的樣機實驗方法，動態(tài)仿真方法具有設(shè)計周期短，設(shè)計成本低等優(yōu)點[1-4]。

制冷系統(tǒng)性能測試需要測試24 h或更長時間，制冷系統(tǒng)實時動態(tài)仿真也相應地需要模擬相同的測試時間下系統(tǒng)的性能，這就需要對蒸發(fā)器模型進行上千萬次的調(diào)用。為了仿真長時間系統(tǒng)性能，蒸發(fā)器模型需要滿足以下要求[2]：

1)穩(wěn)定性：單次計算絕對穩(wěn)定，并保證內(nèi)部質(zhì)量和能量嚴格守恒。否則，一次的計算發(fā)散會導致整個系統(tǒng)仿真的失敗。

2)速度：單次計算時間應該足夠短。如果單次的計算時間過長，則上萬次調(diào)用花費的時間在工程應用中無法接受。

3)精度：系統(tǒng)控制策略需要對蒸發(fā)器的過熱度和冷凝器的過冷度等參數(shù)精確控制，所以要求蒸發(fā)器的模型能夠準確模擬蒸發(fā)器的各狀態(tài)參數(shù)。

現(xiàn)有的蒸發(fā)器動態(tài)仿真模型包括：集中參數(shù)模型、分布參數(shù)模型、移動邊界模型(分相集中參數(shù)模型)。集中參數(shù)模型將蒸發(fā)器看成一個控制單元，具有計算速度快，計算穩(wěn)定等優(yōu)點，但由于采用均相假設(shè)，它不能計算蒸發(fā)器的過熱度等狀態(tài)參數(shù)，計算精度不高[5]。分布參數(shù)模型將蒸發(fā)器分成若干微元計算，具有計算精度高，計算穩(wěn)定等優(yōu)點，但因為涉及大量微元的計算，所以計算速度較慢[6-7]。分布參數(shù)模型的計算速度比移動邊界模型的慢3倍以上[8]。移動邊界模型將蒸發(fā)器的每個相區(qū)作為一個控制單元，它具有計算速度快，精度高，計算穩(wěn)定優(yōu)點，在現(xiàn)有方法中具有綜合優(yōu)勢，它滿足前面列出的蒸發(fā)器動態(tài)仿真所有要求[3-4,9-16].

移動邊界方法雖然是目前計算方法中最可行的方法，但它存在質(zhì)量和能量不嚴格守恒問題[15]，導致采用該方法模擬長時間系統(tǒng)性能時，質(zhì)量和能量的計算誤差在長時間的仿真中誤差不斷積累，影響仿真精度。Cecchhinato L[15]已經(jīng)指出p-ho方法在選取10 s時間步長下，仿真20 h的蒸發(fā)器狀態(tài)，質(zhì)量和能量誤差高達25%以上。

本文的目的是克服現(xiàn)有移動邊界模型質(zhì)量和能量不嚴格守恒問題，提出質(zhì)量和能量嚴格守恒的移動邊界模型。

1 現(xiàn)有移動邊界模型質(zhì)量和能量守恒性分析

1.1 現(xiàn)有蒸發(fā)器移動邊界模型

現(xiàn)有的移動參數(shù)模型根據(jù)選取的標記蒸發(fā)器狀態(tài)的參數(shù)的不同可分為以下三類：

1)p-ho法：選取壓力、出口焓值、兩相區(qū)長度三個參數(shù)作為標記蒸發(fā)器狀態(tài)的參數(shù)，并建立基于選取的狀態(tài)參數(shù)展開質(zhì)量和能量守恒方程，并進行求解[4,10-12]；

三類方法基于各自選取的狀態(tài)參數(shù)展開的質(zhì)量和能量守恒方程，得出各自的控制方程，如表1所示。

表1公式中：p為壓力，Pa;h為出口焓值，kJ/kg;L為換熱管長度，m；A為管橫截面積，m2；γ為空泡系數(shù)；ρ為密度，kg/m3；m為質(zhì)量流量，kg/h；Q為換熱量，W；t為時間，s；Ke為系數(shù)；下標i和o分別表示進口和出口；下標in和out分別表示進入和流出；下標tp和sh分別表示兩相和過熱；下標l和g分別表示飽和液體和飽和氣體，下標int表示相間流動；下標total表示總值；字母上的橫線表示平均值。

以上三類移動邊界模型的狀態(tài)參數(shù)都是壓力項和其他兩項或三項的組合，可統(tǒng)一表示為(N1，N2，N3，N4), 方程可統(tǒng)一表示為公式(1)和(2)，他們都是通過展開能量和質(zhì)量的偏微分項獲得，如公式(3)和(4)所示。

(1)

(2)

(3)

(4)

在公式(1)～公式(4)中，M為蒸發(fā)器總質(zhì)量,kg；E為蒸發(fā)器總能量,kJ；fM為展開的質(zhì)量的偏微分項；fE為展開的能量的偏微分項，它們都是關(guān)于制冷劑的物性的函數(shù);其他參數(shù)的命名與表1公式相同。

1.2 現(xiàn)有模型質(zhì)量和能量的守恒性分析

以上現(xiàn)有控制方程中，質(zhì)量和能量的偏微分項都展開成制冷劑物性的函數(shù)，即公式(1)和(2)中的fM和fE項，而其中的物性函數(shù)必然存在誤差[15]。例如p-ho法質(zhì)量相對壓力的偏微分項中，飽和液體密度和飽和氣體的密度都是通過物性函數(shù)計算，而這些函數(shù)必然存在計算誤差。所以，控制方程中質(zhì)量和能量的偏微分項存在誤差，如公式(5)和公式(6)所示。將公式(5)和公式(6)分別帶入質(zhì)量和能量守恒方程中，可得公式(7)和公式(8)。根據(jù)公式(7)和公式(8)，可導出控制方程質(zhì)量和能量的計算偏差范圍，如公式(9)和公式(10)所示。

(5)

(6)

(7)

(8)

(9)

(10)

由公式(9)和公式(10)可以得出，現(xiàn)有模型控制方程質(zhì)量和能量不嚴格守恒，質(zhì)量和能量的總誤差與仿真時間成正比，與狀態(tài)函數(shù)的步長成正比。產(chǎn)生質(zhì)量和能量不嚴格守恒的根本原因是表示質(zhì)量和能量的偏微分項展開的制冷劑物性的函數(shù)存在誤差。

2 質(zhì)量和能量嚴格守恒移動邊界模型

2.1 基本思路

建立質(zhì)量和能量嚴格守恒的移動邊界模型的基本思路是選取總質(zhì)量M和總能量E作為蒸發(fā)器狀態(tài)參數(shù)，這樣質(zhì)量和能量守恒方程就作為基本控制方程，如公式(11)和公式(12)所示。通過以上方法，控制方程就不會存在誤差，從而質(zhì)量和能量的時間導數(shù)也不會存在誤差，下一時刻的狀態(tài)參數(shù)也不會有誤差，這樣就沒有誤差的傳遞，從而保證了質(zhì)量和能量的嚴格守恒，如圖1所示。

(11)

(12)

通過以上方法可以保證計算中質(zhì)量和能量嚴格守恒，但帶來的問題有：1)如何基于總質(zhì)量和總能量確定制冷劑分布情況；2)如何基于總質(zhì)量和總能量求解制冷劑其他狀態(tài)參數(shù)。

為了解決以上問題，本文提出了：

1)制冷劑不同分布情況下的計算公式；

2)確定制冷劑分布情況以及求解制冷劑其他狀態(tài)參數(shù)的算法。

2.2 計算公式

蒸發(fā)器中制冷劑分布存在兩種狀態(tài)：1)全兩相狀態(tài)；2)兩相和過熱共存狀態(tài)。下面分別列出兩種狀態(tài)下的計算公式。

全兩相情況下，質(zhì)量和能量分別由公式(13)～公式(17)得出。

(13)

(14)

(15)

(16)

(17)

兩相和過熱共存情況下，質(zhì)量和能量分別由公式(18)～公式(21)得出。其中兩相區(qū)平均密度和平均焓值可由公式(15)和公式(17)求出，過熱區(qū)平均密度和平均焓值可由公式(20)和公式(21)求出。

(18)

(19)

(20)

(21)

式(13)～式(21)中:M為蒸發(fā)器的總質(zhì)量，kg;Ecal為蒸發(fā)器的總能量,kJ;x為干度；C為系數(shù)；其他變量的命名與表1相同。

2.3 算法

基于公式(13)～公式(21)，在已知蒸發(fā)器中制冷劑總質(zhì)量和總能量情況下，確定蒸發(fā)器中制冷劑分布狀態(tài)和狀態(tài)參數(shù)的計算方法如下：

1)輸入蒸發(fā)器中制冷劑總質(zhì)量M和總能量E以及蒸發(fā)器進口干度；

2)假設(shè)蒸發(fā)器內(nèi)壓力p；

5)根據(jù)制冷劑不同分布狀態(tài)，分別采用公式(13)～公式(17)或(18)～公式(21)計算制冷劑總能量Ecal；

6)當計算的總能量Ecal和實際總能量E相等時，計算收斂，輸出蒸發(fā)器中制冷劑狀態(tài)參數(shù)，否則重新調(diào)整p，返回步驟3)。

具體算法如圖2所示。

3 模型有效性驗證與對比分析

為了證明新模型有效，本文分別驗證了它的質(zhì)量和能量的守恒性、穩(wěn)定性以及精度和速度。

3.1 質(zhì)量和能量守恒性驗證與對比分析

模型質(zhì)量和能量守恒性的驗證方法是選取典型結(jié)構(gòu)的蒸發(fā)器，輸入的質(zhì)量和能量變化率長時間周期性變化下，考察蒸發(fā)器內(nèi)制冷劑質(zhì)量和能量的守恒性。

用于測試的蒸發(fā)器的結(jié)構(gòu)參數(shù)、初始條件都選擇典型的參數(shù)，制冷劑選取為R134a，仿真步長為1 s，仿真的時間為24 h，如表2所示。輸入的質(zhì)量和能量變化率長時間周期性變化，如圖3和圖4所示。

3.2 穩(wěn)定性驗證與對比分析

本文通過比較48 h后的仿真結(jié)果和初始結(jié)果的偏差，來判別各方法的穩(wěn)定性。圖7～圖11顯示各方法的穩(wěn)定性情況。由圖可以看出，由于輸入的質(zhì)量和能量變化率作周期性變化，蒸發(fā)器制冷劑的狀態(tài)參數(shù)也呈周期性變化，在短時間內(nèi)，各方法的仿真結(jié)果趨于一致，但在仿真到48 h后的性能時，僅本文提出的M-E方法的仿真結(jié)果相對于初始時刻的仿真結(jié)果沒有發(fā)生偏移，而現(xiàn)有其他三種方法的仿真結(jié)果都發(fā)生了偏移，這說明本文提出的M-E方法穩(wěn)定性最好。

3.3 精度和速度驗證與對比分析

采用新方法仿真冰箱蒸發(fā)器正常工作時的性能，仿真結(jié)果和實驗結(jié)果吻合很好，如圖12所示，說明新方法可以用于仿真蒸發(fā)器的動態(tài)性能。

4 結(jié)論

1)對于現(xiàn)有的用于蒸發(fā)器動態(tài)仿真的移動邊界模型，按照制冷劑壓力、兩相區(qū)長度等狀態(tài)參數(shù)展開質(zhì)量和能量方程存在誤差，這是導致了現(xiàn)有的移動邊界模型質(zhì)量和能量不嚴格守恒的根本原因。

2)選取制冷劑總質(zhì)量和總能量作為狀態(tài)參數(shù)，將質(zhì)量和能量守恒方程直接作為控制方程，能夠避免方程展開產(chǎn)生的誤差，從而保證了新方法的質(zhì)量和能量嚴格守恒。

3)為了實現(xiàn)在已知總質(zhì)量和總能量的情況下對蒸發(fā)器中制冷劑分布情況和其他狀態(tài)參數(shù)的求解，本文給出蒸發(fā)器制冷劑在不同分布情況下的計算公式和計算方法。

4)通過仿真案例可以得出，在質(zhì)量和能量守恒方面，新方法質(zhì)量和能量嚴格守恒，而現(xiàn)有其他方法都不能同時保證質(zhì)量和能量的嚴格守恒；在穩(wěn)定性方面，新方法的仿真結(jié)果沒有發(fā)生任何偏移，具有良好的穩(wěn)定性，而其他方法的仿真結(jié)果都發(fā)生了偏移現(xiàn)象；在計算精度和速度方面，新方法計算精度和已有方法相當，新方法計算速度略低。

[1] 邵雙全, 石文星, 陳華俊, 等. 制冷空調(diào)系統(tǒng)計算機仿真技術(shù)綜述[J]. 制冷與空調(diào), 2002, 2(3): 10-15. (Shao Shuangquan, Shi Wenxing, Chen Huajun, et al. Literature review of computer simulation on refrigeration and air conditioning system[J]. Refrigeration and Air-Conditioning, 2002, 2(3): 10-15.)

[2] Ding Guoliang. Recent developments in simulation techniques for vapour-compression refrigeration systems[J]. International Journal of Refrigeration, 2007, 30(7): 1119-1133.

[3] Leducq D, Guilpart J, Trystram G. Non-linear predictive control of a vapour compression cycle[J]. International Journal of Refrigeration, 2006, 29(5): 761-772.

[4] Rasmussen B P, Alleyne A G. Control-oriented modeling of transcritical vapor compression systems[J]. Journal of Dynamic Systems, Measurement and Control, 2004, 126: 54-64.

[5] Chi J, Didion D. A simulation model of the transient performance of a heat pump[J]. International Journal of Refrigeration, 1982, 5(3): 176-184.

[6] Gruhle W D, Isermannn R. Modeling and control of a refrigerant evaporator[J]. Journal of Dynamic Systems, Measurement and Control, 1985, 107: 235-240.

[7] MacArthur J W, Grald E W. Unsteady compressible two-phase flow model for predicting cyclic heat pump performance and a comparison with experimental data[J]. International Journal of Refrigeration, 1989, 12(1): 29-41.

[8] Bendapudi S, Braun J E, Groll E A. A comparison of moving-boundary and finite-volume formulations for transients in centrifugal chillers[J]. International Journal of Refrigeration, 2008, 31(8): 1437-1452.

[9] He X D, Liu S, Asada H. Modeling of vapor compression cycles for multivariable feedback control of HVAC systems[J]. Journal of Dynamic Systems, Measurement and Control, 1997, 119(2): 183-191.

[10] Pettit N B O L, Willatzen M, Ploug-Sorensen L. A general dynamic simulation model for evaporators and condensers in refrigeration part two: simulation and control of an evaporator[J]. International Journal of Refrigeration, 1998, 21(5): 404-414.

[11] Zhang Weijiang, Zhang Chunlu. A generalized moving-boundary model for transient simulation of dry-expansion evaporators under larger disturbances[J]. International Journal of Refrigeration, 2006, 29(7): 1119-1127.

[12] Liang N, Shao S, Tian C, et al. Dynamic simulation of variable capacity refrigeration systems under abnormal conditions[J]. Applied Thermal Engineering, 2010, 30(10): 1205-1214.

[13] Li B, Alleyne A G. A dynamic model of a vapor compression cycle with shut-down and start-up operations[J]. International Journal of Refrigeration, 2010, 33(3): 538-552.

[14] McKinley T L, Alleyne A G. An advanced nonlinear switched heat exchanger model for vapor compression cycles using the moving-boundary method[J]. International Journal of Refrigeration, 2008, 31(7): 1253-1264.

[15] Cecchinato L, Mancini F. An intrinsically mass conservative switched evaporator model adopting the moving-boundary method[J].International Journal of Refrigeration, 2012, 35(2): 349-364.

[16] Lei Z, Zaheeruddin M. Dynamic simulation and analysis of a water chiller refrigeration system[J]. Applied Thermal Engineering, 2005, 25(14/15): 2258-2271.

About the corresponding author

Ding Guoliang, male, Ph. D. / professor, Institute of Refrigeration and Cryogenics, Shanghai Jiaotong University, +86 21-34203278, E-mail: glding@sjtu.edu.cn. Research fields: simulation and optimization research for room air conditioner and utilization of new refrigerant.

An Intrinsically Mass and Energy Conservative Model for Transient Simulation of Evaporator

Zhao Dan Ding Guoliang Hu Haitao

(Institute of Refrigeration and Cryogenics, Shanghai Jiaotong University, Shanghai, 200240, China)

The existing moving-boundary models for evaporators are obtained by expanding the continuity and internal energy equations, and the error caused by the expansion results that the models are not intrinsically conservative on both mass and energy, resulting accuracy deterioration. To guarantee intrinsically conservative on both mass and energy in calculation, the mass and internal energy are selected as the state variables, and the continuity and internal energy equations are directly solved to avoid the expansion error. The equations and algorithms for refrigerant state and phase distribution in evaporator are developed when the mass and internal energy are given. The case study shows the new model is steady and intrinsically conservative on both mass and energy in the performance simulation of 48 hours, and the simulation results agree well with the experimental data.

evaporator; transient simulation; moving-boundary model

0253- 4339(2015) 01- 0076- 08

10.3969/j.issn.0253- 4339.2015.01.076

2014年4月20日

TB61+5； TB657.5； TP391.9

A

丁國良，男，教授，博導，上海交通大學機械與動力工程學院制冷所，(021) 34206378，E-mail：glding@sjtu.edu.cn。研究方向：制冷空調(diào)裝置的仿真、優(yōu)化與新工質(zhì)應用。