段鐘弟，任滔，丁國良，陳杰，浦暉，密曉光

（1上海交通大學(xué)制冷與低溫工程研究所，上海200240；2中海石油氣電集團(tuán)技術(shù)研發(fā)中心，北京100028）

引 言

天然氣作為一種高效清潔的能源，近年來得到了廣泛的應(yīng)用［1］。我國南海擁有豐富的待開發(fā)天然氣資源，已探明的儲量占我國天然氣總儲量的1／3［2］。開采海上天然氣的基礎(chǔ)裝置是浮式天然氣生產(chǎn)儲卸平臺 （LNG－FPSO）。LNG－FPSO的主低溫?fù)Q熱器為繞管式換熱器，具有占地面積小、力學(xué)性能好和抗泄漏能力強(qiáng)等優(yōu)點(diǎn)［3］。

繞管式換熱器由數(shù)以千計(jì)的管子按層排布螺旋纏繞在中心筒上，單個(gè)換熱器投資巨大，不能基于實(shí)驗(yàn)進(jìn)行設(shè)計(jì)，需要先進(jìn)的設(shè)計(jì)方法［4］?；诜抡娴脑O(shè)計(jì)方法具有高效、低成本的優(yōu)勢，能夠應(yīng)用于繞管式換熱器設(shè)計(jì)［5］。LNG－FPSO是海上浮動平臺，會隨風(fēng)浪在海上晃動，從而造成繞管式換熱器晃動，使得其內(nèi)部流體的溫度、壓力隨時(shí)間發(fā)生波動［6］，因此需要能夠預(yù)測繞管式換熱器動態(tài)特性的模型。

因天然氣液化獨(dú)特的工藝和需求，繞管式換熱器的結(jié)構(gòu)和工作具有如下特征，包括：①工質(zhì)是由烴類和氮?dú)饨M成的混合冷劑；②工況范圍覆蓋包含亞臨界和超臨界的所有相區(qū)；③多股流體并行換熱；④單個(gè)換熱器內(nèi)冷熱流體并發(fā)相變，同時(shí)存在蒸發(fā)和冷凝?；诩倕?shù)的多股流換熱器動態(tài)模型［7－10］不能反映繞管式換熱器中的相變特征；基于分布參數(shù)的多股流動態(tài)模型［11－12］涉及大量微元的計(jì)算，計(jì)算耗時(shí)長［13］。因此，本研究采用具有良好精度和計(jì)算速度的移動邊界模型［14］來預(yù)測繞管式換熱器的動態(tài)特性。現(xiàn)有的移動邊界模型可以對繞管式換熱器前兩項(xiàng)特征中的多組分工質(zhì)、多相區(qū)模式進(jìn)行計(jì)算，具體見表1；而對于后兩項(xiàng)特征中多股流并行換熱和冷熱流體并發(fā)相變的計(jì)算，還需要進(jìn)一步的研究。

表1 已有研究中代表性的換熱器動態(tài)模型Table 1 Existing representative dynamic models of heat exchangers

多股流并行換熱中，一股流體同時(shí)與其他多股流體發(fā)生熱交換，使得各股冷熱流體之間的換熱關(guān)系復(fù)雜化［21］。同時(shí)，在移動邊界模型下流體的控制單元按照相區(qū)進(jìn)行劃分，由于各個(gè)相區(qū)的長度不一，在空間分布上互相交錯，使得分屬不同流路的相區(qū)之間換熱關(guān)系復(fù)雜化。以上兩方面的原因?qū)е乱苿舆吔缒Ｐ椭械亩喙闪鞑⑿袚Q熱難以計(jì)算。針對特定多股流形式的建模方法［8－9］不能滿足繞管式換熱器各種不同流程布置的要求；采用基于集總參數(shù)的互連矩陣來描述流路換熱關(guān)系的方法［7］不適用于移動邊界模型下相區(qū)之間的換熱計(jì)算。因此，需要一種能夠描述繞管式換熱器多股流多相區(qū)的換熱關(guān)系，對任意的流路和相區(qū)分布形式進(jìn)行計(jì)算的方法。

當(dāng)換熱器中同時(shí)存在蒸發(fā)和冷凝相變過程時(shí)，管壁兩側(cè)的冷熱流體均存在相邊界。在動態(tài)過程中，兩側(cè)流體的相邊界同時(shí)發(fā)生移動，一方面會出現(xiàn)相邊界的 “消失”和 “重現(xiàn)”，另一方面還存在冷熱流體相邊界的 “交錯”，使得邊界的移動行為更加復(fù)雜。相邊界的 “消失”和 “重現(xiàn)”要求移動邊界模型能夠處理不同的相區(qū)模式并能夠光滑轉(zhuǎn)換，這一部分已有大量的文獻(xiàn)進(jìn)行了研究［15－20，22－29］；相邊界的 “交錯”則會使得相邊界的排列形式多種多樣。在不同的排列形式下，管壁兩側(cè)的冷熱流體相區(qū)分布會發(fā)生變化，使得傳熱計(jì)算的邊界條件不同。因此，需要一種能夠計(jì)算熱交換兩側(cè)流體均存在相邊界時(shí)的傳熱，并且能夠處理相邊界 “交錯”的情形的方法。

因此，要建立繞管式換熱器移動邊界模型，需要：①多股流多相區(qū)的換熱關(guān)系的數(shù)學(xué)描述方法；②建立相邊界交錯時(shí)的換熱計(jì)算模型。

1 多股流多相區(qū)換熱關(guān)系的

LNG繞管式換熱器通常為級聯(lián)形式，級數(shù)目根據(jù)液化流程有所區(qū)別，換熱器結(jié)構(gòu)如圖1所示。

圖1 繞管式換熱器級聯(lián)示意圖［5］Fig.1 Schematic diagram of SWHE bundles［5］

在單級管束中，管側(cè)流體自管束底部向上流動，殼側(cè)流體自管束頂部向下流動。管束內(nèi)部各個(gè)相區(qū)長度不一，單個(gè)相區(qū)同時(shí)與多個(gè)相區(qū)發(fā)生換熱，如圖2所示。

圖2 單級管束相區(qū)分布示意圖Fig.2 Distribution of phase regions in single bundle

為了描述一個(gè)具有n股管側(cè)流路、1股殼側(cè)流體且每股流路中有多個(gè)相區(qū)的換熱器，可以通過二維矩陣的形式表示換熱計(jì)算中的相關(guān)參數(shù)。二維矩陣的行和列數(shù)目相同，一共由3（n＋1）×3（n＋1）個(gè)元素組成，每一行 （列）代表一股流路的單個(gè)相區(qū)。矩陣中的元素上角標(biāo)表示當(dāng)前行對應(yīng)的流路編號和相區(qū)編號，下角標(biāo)表示當(dāng)前列對應(yīng)的流路編號和相區(qū)編號。對于繞管式換熱器，殼程流路對應(yīng)的編號為＃0，管程流路數(shù)目為n時(shí)對應(yīng)的流路編號依次為＃1～＃n。每股流路中統(tǒng)一將過熱區(qū)編號為＃1，兩相區(qū)編號為＃2，過冷區(qū)編號為＃3。天然氣流路的超臨界區(qū)視為過熱區(qū)，編號為＃1。二維矩陣的表達(dá)形式如圖3所示。

圖3 多股流多相區(qū)二維矩陣Fig.3 Matrix of multiple streams and multiple zones

該二維矩陣可以視作一個(gè)分塊矩陣，由 （n＋1）×（n＋1）個(gè)子矩陣構(gòu)成，每一個(gè)子矩陣由3×3個(gè)元素構(gòu)成，子矩陣的表達(dá)形式如式 （1）所示。

矩陣中的元素可以用換熱計(jì)算中具體的物理參數(shù)代替。這些參數(shù)可以分為兩類：一類是與相區(qū)分布特征無關(guān)的參數(shù)，包括流體溫度Tf、單位長度換熱面積Funit、傳熱系數(shù)；另一類是與相區(qū)分布特征相關(guān)的參數(shù)，包括管壁溫度Tw、換熱長度L以及換熱量Q。

當(dāng)矩陣表示與相區(qū)分布特征無關(guān)的參數(shù)時(shí)，矩陣為對角矩陣，只有在矩陣對角線上的元素具有物理意義，代表當(dāng)前流路當(dāng)前相區(qū)的物理參數(shù)。例如，當(dāng)矩陣表示流體溫度時(shí)，矩陣表達(dá)形式如式（2）所示。

其中，子矩陣的表達(dá)形式如式 （3）所示。

當(dāng)矩陣表示與相區(qū)分布特征相關(guān)的參數(shù)時(shí)，矩陣為互連矩陣，只有當(dāng)流路i和流路j存在換熱關(guān)系時(shí)對應(yīng)的子矩陣中的元素才有值。對角線上的子矩陣以及其他行列對應(yīng)的流路間不存在換熱關(guān)系的子矩陣，里面的元素將沒有值。例如，當(dāng)矩陣表示換熱長度時(shí)，矩陣表達(dá)形式如式 （4）所示。

其中，子矩陣Lij的表達(dá)形式如式 （5）所示。

這里使用的長度l為歸一化長度，表示當(dāng)前段長度占總長度的比值。在長度矩陣中，任一元素）具有的物理意義為：流路i的第m個(gè)相區(qū)與流路j的第n個(gè)相區(qū)之間的換熱長度。通過單個(gè)長度子矩陣可以描述兩股流路中的相區(qū)分布情況，如圖4所示。

圖4 管殼側(cè)兩股流路的相區(qū)分布示意圖Fig.4 Distribution of phase regions in tube and shell side

對應(yīng)圖4中的相區(qū)分布情況，長度子矩陣的表達(dá)形式如式 （6）所示。

當(dāng)矩陣中的元素表示管壁溫度時(shí)，其物理意義為：流路i的第m個(gè)相區(qū)與流路j的第n個(gè)相區(qū)之間的管壁區(qū)間的平均溫度。對于圖3的相區(qū)分布情況，溫度子矩陣Tiw0的表達(dá)形式如式 （7）所示。

2 相邊界交錯下的傳熱計(jì)算模型

相變界交錯下的傳熱計(jì)算模型包括流體的控制方程、流體與管壁的傳熱方程和管壁的控制方程。其中，流體的控制方程根據(jù)相區(qū)的質(zhì)量和能量守恒建立，采用 的 移 動 邊 界 模 型 的 假 設(shè)［15，17－18，20］包 括：①流體的流動視為一維流動；②單股流體內(nèi)壓力均勻分布，忽略流動過程中的壓力損失；③忽略流體的軸向?qū)?；④忽略管壁?nèi)部的導(dǎo)熱。

根據(jù)上述假設(shè)，微分形式的流體控制方程如式（8）和式 （9）所示。

式中，ρ表示密度，A表示流體的流通面積，m表示質(zhì)量流量，h表示比焓，p表示壓力，Q表示流體的換熱量，τ表示時(shí)間，z表示流體的流動方向。

將上述偏微分方程沿流體的相區(qū)長度積分，可以得到不同展開形式的微分代數(shù)方程［30］。對展開后的微分代數(shù)方程進(jìn)行求解，可以得到流體相區(qū)的壓力、溫度等狀態(tài)參數(shù)以及相區(qū)長度隨時(shí)間的變化。對于相變界交錯下的多股流繞管式換熱器，求解上式的關(guān)鍵在于建立單個(gè)流體相區(qū)與管壁之間的傳熱模型和管壁溫度的計(jì)算模型。

2.1 流體與管壁的傳熱模型

流體與管壁的傳熱模型，可以根據(jù)第1節(jié)中建立的流體溫度、管壁溫度、換熱長度的二維矩陣，使用統(tǒng)一的表達(dá)形式進(jìn)行計(jì)算。

2.1.1 兩個(gè)相區(qū)之間的換熱計(jì)算 假定相區(qū)A為流路i的第m個(gè)相區(qū)，相區(qū)B為流路j的第n個(gè)相區(qū)。由相區(qū)A傳遞給相區(qū)B的換熱量可以通過式（10）進(jìn)行計(jì)算

反之，由相區(qū)B傳遞給相區(qū)A的熱量為

2.1.2 單個(gè)相區(qū)的總換熱量計(jì)算 假定相區(qū)A為流路i的第m個(gè)相區(qū)，則相區(qū)A的總換熱量為相區(qū)A與其他所有相區(qū)的換熱量之和，可以通過式（12）進(jìn)行計(jì)算

2.2 管壁溫度的計(jì)算模型

熱交換兩側(cè)的冷熱流體均發(fā)生相變時(shí)，管壁邊界的排列有多種形式，圖5是其中兩種代表性的邊界排布示意圖。對于管壁邊界，最多可能出現(xiàn)的邊界數(shù)為6個(gè)，如圖5（a）所示。對于只有兩個(gè)相區(qū)或者單個(gè)相區(qū)的情況，保留所有的邊界編號，對于其中實(shí)際不存在的相區(qū)，將其邊界視為重疊在一起，如圖5（b）所示。

圖5 管壁邊界排布示意圖Fig.5 Boundary arrangement of tube wall

對于其中的單個(gè)管壁區(qū)間 （區(qū)間k），其能量傳遞過程如圖6所示。根據(jù)圖6，并應(yīng)用萊布尼茲積分法則，得到微分形式的控制方程，如式 （13）所示。

圖6 管壁熱量傳遞示意圖Fig.6 Heat transfer in tube wall

式 （19）中，lw，m＋n、Tw，m＋n根據(jù)式 （20）計(jì)算得到

為了使模型能夠在相邊界發(fā)生 “交錯”時(shí)光滑切換，虛擬區(qū)間的管壁溫度應(yīng)保持與對應(yīng)實(shí)際區(qū)間的管壁溫度一致，通過式 （21）進(jìn)行計(jì)算。

式中，Kw為松弛因子，Tw，track為對應(yīng)的實(shí)際區(qū)間的管壁溫度值。

3 計(jì)算實(shí)例和驗(yàn)證

3.1 計(jì)算實(shí)例

選取一個(gè)典型流路布置的繞管式換熱器進(jìn)行實(shí)例計(jì)算，其流路方案如圖7所示。

設(shè)置具體的工藝方案：流路1和流路2為熱流體，流路3為冷流體；流路1的天然氣運(yùn)行在超臨界區(qū)，流路2的混合冷劑運(yùn)行在亞臨界區(qū)并由氣態(tài)轉(zhuǎn)變?yōu)閮上?，流?的殼側(cè)冷劑運(yùn)行在亞臨界區(qū)并由兩相轉(zhuǎn)變?yōu)闅鈶B(tài)。初始穩(wěn)態(tài)時(shí)的具體工藝方案見表2。

圖7 計(jì)算實(shí)例的流路方案Fig.7 Flow arrangement in case study

表2 計(jì)算實(shí)例的初始工藝方案Table 2 Initial parameters of case

設(shè)置擾動條件為：殼側(cè)流體的傳熱系數(shù)隨時(shí)間周期性波動。實(shí)例計(jì)算中流體與管壁的傳熱計(jì)算、管壁溫度的計(jì)算采用本研究提出的模型，流體相區(qū)的控制方程采用文獻(xiàn)中的模型［19］，傳熱系數(shù)采用換熱關(guān)聯(lián)式進(jìn)行計(jì)算 （表3）。

表3 計(jì)算實(shí)例中的換熱關(guān)聯(lián)式Table 3 List of applied heat transfer correlations

3股流體的壓力、出口溫度和相區(qū)長度隨時(shí)間的變化分別如圖8～圖10所示。

圖8 各股流路的壓力隨時(shí)間的變化曲線Fig.8 Graph of pressure in each stream with time

圖9 各股流路的出口溫度隨時(shí)間的變化曲線Fig.9 Graph of outlet temperature in each stream with time

圖10 各股流路的相區(qū)長度隨時(shí)間的變化曲線Fig.10 Graph of zone length in each stream with time

3股流體中，流路2和流路3均發(fā)生相變，其相變界面處的邊界隨時(shí)間移動的曲線如圖11所示。圖中，兩條曲線發(fā)生交叉時(shí)，表示流路2和流路3的相邊界發(fā)生 “交錯”的情況。

圖11 相界面隨時(shí)間的移動曲線Fig.11 Graph of interface move with time

3.2 精度驗(yàn)證

目前公開文獻(xiàn)中尚未發(fā)表過LNG或FLNG繞管式換熱器的相關(guān)實(shí)驗(yàn)數(shù)據(jù)和實(shí)際工廠的動態(tài)運(yùn)行數(shù)據(jù)，因此沒有實(shí)驗(yàn)數(shù)據(jù)可以用來與模型比對?？紤]到本研究做的創(chuàng)新工作主要集中在如何描述和計(jì)算多股流的多相區(qū)耦合的問題，對于其中的任意一股流體，本研究用的模型與已有文獻(xiàn)上的分相移動邊界模型一致，因此本研究開發(fā)的模型的正確性可以通過已有的單股流的分相移動邊界模型來驗(yàn)證。

驗(yàn)證時(shí)流體建模采用與文獻(xiàn)中模型相同的P－ˉγ法，流體與管壁的換熱計(jì)算和管壁溫度的計(jì)算采用本研究提出的模型。驗(yàn)證實(shí)例采用如下工藝方案：管側(cè)為單股流體形式，工質(zhì)采用C2和C3組成的混合冷劑，并設(shè)置管側(cè)為過熱氣體、殼側(cè)為兩相到過熱的相變流體。驗(yàn)證時(shí)傳熱系數(shù)取為常數(shù)，這樣避免了關(guān)聯(lián)式誤差對驗(yàn)證精度的影響，但在實(shí)際應(yīng)用中可以采用表3列出的換熱關(guān)聯(lián)式或其他關(guān)聯(lián)式進(jìn)行計(jì)算。空泡系數(shù)模型采用均相模型，同樣可以由其他模型代替。流體的物性計(jì)算采用GERG2008［35］方 法，通 過 REFPROP 6.01 生 成。具體的工藝參數(shù)和結(jié)構(gòu)參數(shù)見表4。

設(shè)置擾動條件為：在初始時(shí)換熱器處于穩(wěn)態(tài)，隨后保持管側(cè)和殼側(cè)的進(jìn)口焓值不變，增大殼側(cè)流量，此處增量取為5%。

仿真結(jié)果顯示，兩種模型下殼側(cè)及管側(cè)的壓力和出口焓值的計(jì)算結(jié)果吻合良好，最大偏差小于4%，如圖12和圖13所示。

表4 繞管式換熱器精度驗(yàn)證參數(shù)Table 4 SWHE parameters of validation case

圖12 殼側(cè)壓力和出口焓值隨時(shí)間的變化曲線Fig.12 Graph of shell－side pressure and outlet enthalpy with time

圖13 管側(cè)壓力和出口焓值隨時(shí)間的變化曲線Fig.13 Graph of tube－side pressure and outlet enthalpy with time

4 結(jié) 論

（1）建立的二維矩陣數(shù)學(xué)描述方法能夠?qū)@管式換熱器多股流路中各個(gè)相區(qū)之間的換熱關(guān)系進(jìn)行描述，實(shí)現(xiàn)了移動邊界下多股流并行換熱的計(jì)算。

（2）建立的傳熱計(jì)算模型能夠?qū)錈崃黧w并發(fā)相變的情況進(jìn)行計(jì)算，并可以處理相邊界的 “交錯”，實(shí)現(xiàn)了移動邊界下蒸發(fā)和冷凝并存的計(jì)算。

（3）實(shí)例表明本研究提出的建模方法可以應(yīng)用于多股流繞管式換熱器的動態(tài)仿真，與文獻(xiàn)中模型的仿真結(jié)果吻合良好，最大偏差小于4%，滿足精度要求。

符 號 說 明

Aw——管壁橫截面積，m2

Cw——管壁比熱容，J·kg－1·K－1

Funit——單位長度換熱面積，m2·m－1

Kw——時(shí)間常數(shù)的倒數(shù)，s－1

L——換熱長度矩陣

L——總換熱長度，m

l——?dú)w一化的換熱長度

lw——?dú)w一化的管壁區(qū)間長度

Q——換熱量，W

Tf——流體溫度矩陣

Tf——流體溫度，K

Tint——管壁區(qū)間交界面溫度，K

Tw——管壁溫度矩陣

Tw——管壁溫度，K

α——傳熱系數(shù)，W·m－2·K－1

ρw——管壁材料密度，kg·m－3

τ——時(shí)間，s

χl——管壁區(qū)間的下邊界位置

χu——管壁區(qū)間的上邊界位置

χw——管壁區(qū)間的邊界位置

下角標(biāo)

f——流體

in——向內(nèi)傳遞的

int——管壁區(qū)間交界面

i，j——流路編號

k——管壁區(qū)間編號

l——管壁區(qū)間的下邊界

m，n——相區(qū)編號

out——向外傳遞的

unit——單位長度

w——管壁

［1］ Momirlan M，Veziroglu T N.The properties of hydrogen as fuel tomorrow in sustainable energy system for a cleaner planet ［J］.InternationalJournalofHydrogenEnergy，2005，30 （7）：795－802.

［2］ Zhang Fengjiu （張鳳久）.Outlook for natural gas development in the South China Sea ［J］.Natur.GasInd.（天然氣工業(yè)），2009，29 （1）：17－20.

［3］ Bukowski J，Liu Y N，Boccella S，etal.Innovations in natural gas liquefaction technology for future LNG plants and floating LNG facilities／／Proceedings of the International Gas Union Research Conference ［C］.2011.

［4］ Pacio J C，Dorao C A.A review on heat exchanger thermal hydraulic models for cryogenic applications ［J］.Cryogenics，2011，51 （7）：366－379.

［5］ Wang T，Ding G，Duan Z，etal.A distributed－parameter model for LNG spiral wound heat exchanger based on graph theory ［J］.AppliedThermalEngineering，2015，81（25）：102－113.

［6］ Lex T，Ohlig K，fredheim A O，etal.Ready for floating LNG－qualification of spiral wound heat exchangers／／Proceedings 15th International Conferences on Liquefied Natural Gas ［C］.Barcelona，Spain，2007：PO－27.1－PO－27.14.

［7］ Luo X，Guan X，Li M，etal.Dynamic behaviour of onedimensional flow multistream heat exchangers and their networks ［J］.InternationalJournalofHeatandMass Transfer，2003，46 （4）：705－715.

［8］ Mishra M，Das P K，Sarangi S.Dynamic behavior of threefluid crossflow heat exchangers ［J］.JournalofHeat Transfer，2008，130 （1）：118－125.

［9］ Das S K，Murugesan K.Transient response of multipass plate heat exchangers with axial thermal dispersion in fluid［J］.InternationalJournalofHeatandMassTransfer，2000，43 （23）：4327－4345.

［10］ Zhang Qin （張勤），Li Zhixin （李志信），Liang Xingang（梁新剛）.Field synergy analysis of the dynamic process of a multi－stream heat exchanger ［J］.JournalofEngineering forThermalEnergy＆Power（熱能動力工程），2009，24（6）：782－786.

［11］ Hammer M.Dynamic simulation of a natural gas liquefaction plant ［D］.Trondheim：Norwegian University of Science and Technology，2004：55－108.

［12］ Xu Yifeng （徐益峰），Cai Zuhui（蔡祖恢）.Research on dynamic model of multi－fluid plate－fin heat exchangers ［J］.JournalofChemicalIndustryandEngineering（China）（化工學(xué)報(bào)），1998，49 （6）：721－728.

［13］ Lei Z，Zaheeruddin M.Dynamic simulation and analysis of a water chiller refrigeration system ［J］.AppliedThermal Engineering，2005，25 （14）：2258－2271.

［14］ Ding G L.Recent developments in simulation techniques for vapour－compression refrigeration systems ［J］.International JournalofRefrigeration，2007，30 （7）：1119－1133.

［15］ Munch Jensen J，Tummescheit H.Moving boundary models for dynamic simulations of two－phase flows／／Modelica’2002 Proceedings ［C］.2002：235－244.

［16］ Rasmussen B P，Alleyne A G.Control－oriented modeling of transcritical vapor compression systems ［J］.Journalof DynamicSystems，Measurement，andControl，2004，126（1）：54－64.

［17］ Zhang W J，Zhang C L.A generalized moving－boundary model for transient simulation of dry－expansion evaporators under larger disturbances ［J］.InternationalJournalof Refrigeration，2006，29 （7）：1119－1127.

［18］ McKinley T L，Alleyne A G.An advanced nonlinear switched heat exchanger model for vapor compression cycles using the moving－boundary method ［J］.InternationalJournalof Refrigeration，2008，31 （7）：1253－1264.

［19］ Li B，Alleyne A G.A dynamic model of a vapor compression cycle with shut－down and start－up operations ［J］.InternationalJournalofRefrigeration，2010，33 （3）：538－552.

［20］ Cecchinato L，Mancini F.An intrinsically mass conservative switched evaporator model adopting the moving－boundary method ［J］.InternationalJournalofRefrigeration，2012，35 （2）：349－364.

［21］ Picon－Nunez M，Polley G T，Medina－Flores M.Thermal design of multi－stream heat exchangers ［J］.Applied ThermalEngineering，2002，22 （14）：1643－1660.

［22］ Willatzen M，Pettit N，Ploug－S?rensen L.A general dynamic simulation model for evaporators and condensers in refrigeration （Ⅰ ）： Moving－boundary formulation of twophase flows with heat exchange ［J］.InternationalJournal ofRefrigeration，1998，21 （5）：398－403.

［23］ Ding G，Zhang C，Lu Z.Dynamic simulation of natural convection bypass two－circuit cycle refrigerator－freezer and its application（Ⅰ）：Component models ［J］.AppliedThermalEngineering，2004，24 （10）：1513－1524.

［24］ Fu L，Ding G，Zhang C.Dynamic simulation of air－to－water dual－mode heat pump with screw compressor ［J］.Applied ThermalEngineering，2003，23 （13）：1629－1645.

［25］ Schurt L C，Hermes C J L，Neto A T.A model－driven multivariable controller for vapor compression refrigeration systems ［J］.InternationalJournalofRefrigeration，2009，32 （7）：1672－1682.

［26］ Pettit N，Willatzen M，Ploug－S?rensen L.A general dynamic simulation model for evaporators and condensers in refrigeration （Ⅱ）：Simulation and control of an eva－porator［J］.InternationalJournalofRefrigeration，1998，21 （5）：404－414.

［27］ Beghi A，Cecchinato L.A simulation environment for dryexpansion evaporators with application to the design of autotuning control algorithms for electronic expansion valves［J］.InternationalJournalofRefrigeration，2009，32 （7）：1765－1775.

［28］ Liang N，Shao S，Tian C，etal.Dynamic simulation of variable capacity refrigeration systems under abnormal conditions ［J］.AppliedThermalEngineering，2010，30（10）：1205－1214.

［29］ Bendapudi S，Braun J E，Groll E A.A comparison of moving－boundary and finite－volume formulations for transients in centrifugal chillers ［J］.InternationalJournal ofRefrigeration，2008，31 （8）：1437－1452.

［30］ Zhao Dan （趙 丹），Ding Guoliang （丁 國 良），Hu Haitao（胡 海 濤 ）.An intrinsically mass and energy conservative model for transient simulation of evaporator ［J］.Journalof Refrigeration（制冷學(xué)報(bào)），2015，36 （1）：76－83.

［31］ Shah M.A general correlation for heat transfer during film condensation inside pipes ［J］.InternationalJournalofHeat＆MassTransfer，1979，22 （4）：547－556.

［32］ Boyko L D，Kruzhilin G N.Heat transfer and hydraulic resistance during condensation of steam in a horizontal tube and in a bundle of tubes ［J］.InternationalJournalofHeat＆MassTransfer，1967，10 （67）：361－373.

［33］ Neeraas B O，F(xiàn)redheim A O，Aunan B.Experimental shellside heat transfer and pressure drop in gas flow for spiralwound LNG heat exchanger ［J］.InternationalJournalof Heat＆MassTransfer，2004，47 （2）：353－361.

［34］ Neeraas B O，F(xiàn)redheim A O，Aunan B.Experimental data and model for heat transfer，in liquid falling film flow on shell－side，for spiral－wound LNG heat exchanger ［J］.InternationalJournalofHeat＆MassTransfer，2004，47（14）：3565－3572.

［35］ Kunz O，Wagner W.The GERG－2008wide－range equation of state for natural gases and other mixtures：an expansion of GERG－2004 ［J］.JournalofChemical＆Engineering Data，2012，57 （11）：3032－3091.