基于一維分布參數法的螺旋管蒸發器動態建模

閔智勇, 黨建軍, 蘇 浩, 張佳楠, 李代金

?

基于一維分布參數法的螺旋管蒸發器動態建模

閔智勇, 黨建軍, 蘇 浩, 張佳楠, 李代金

(西北工業大學航海學院, 陜西西安, 710072)

金屬燃料能源因為其較高的能量密度成為當今水下熱動力系統的一個重要發展方向?；诖? 文中針對Li/SF6熱管反應器內部的螺旋管蒸發器, 采用一維分布參數法建立動態模型, 對不同輸入參數變化進行了仿真計算, 并將結果與移動邊界法計算結果進行對比。分析表明, 一維分布參數模型所得結果相比移動邊界模型與實際系統具有更高的貼合度, 能夠更為準確地反映系統的動態響應過程, 并通過合適的公式推導和算法設計, 得到與移動邊界模型接近的計算效率, 在輸入參數劇烈變化時, 依舊顯示出了良好的模型強壯性。文中研究內容可為水下熱動力系統設計、控制研究和實時系統仿真提供依據。

水下熱動力系統; 螺旋管蒸發器; 一維分布參數法; 移動邊界法; 動態建模

0 引言

美軍MK50魚雷以高能量密度的Li/SF6為熱源[1], 采用朗肯循環構建了一套水下閉式循環熱動力系統, 堪稱典范。此后, Blakeslee[2]等人針對Li/SF6在金屬絲網表面的緩慢燃燒展開研究, 提出了Li/SF6熱管反應器的概念, 旨在將該能源應用于UUV動力系統之上。該動力系統基本構型為Li/SF6在密閉的熱管反應器中反應釋熱, 使液體鋰蒸發, 安裝在反應器外殼中的螺旋管蒸發器吸收Li蒸汽的凝結潛熱, 將工質水加熱至過熱狀態推動渦輪機做功, 乏汽經回熱器加熱進入螺旋管蒸發器的工質水, 再進入冷凝器冷凝后進入集液器, 完成循環過程。

目前, 國內對于Li/SF6熱管反應器的研究主要集中在燃燒反應過程仿真[3-4]和試驗研究[5]等方面, 而對其內部螺旋管蒸發器的研究則顯得相對不足。但在該熱動力系統中, 螺旋管蒸發器將工質水加熱為過熱蒸汽, 其性能優良直接決定渦輪機做功能力, 從而影響整體系統性能。故針對熱管反應器內部存在相變的螺旋管蒸發器, 建立精確的數學模型, 展現其變工況工作時各參數動態響應過程, 對該動力系統的設計及控制有著重要意義。

對于存在相變的管式蒸發器動態建模, 普遍存在2種方法: 一種是移動邊界法, 此方法早期由He X[6]提出, 該方法的核心思想是通過對不同相區分別建立集中參數模型, 以得到優于一般集中參數模型精度的優點, 同時保證了較高的計算效率, 該方法優化、改進后被應用于各相關領域[7-8]。白杰[9]借鑒該方法對Li/SF6熱管反應器內部的螺旋管蒸發器建立了動態模型, 得到了入口流量變化后的響應規律, 但對于結果出現的反沖現象卻無法解釋, 并且壁面條件過于簡化忽略了管壁的熱容效應; 另一種方法是一維分布參數法[10-11], 該方法將蒸發器劃分為多個控制體, 直接對偏微分方程進行離散, 故可反應工質沿流動方向的參數分布變化規律, 具有極高的精度。但限于計算效率較低, 不適用于系統級仿真而常被忽視。

文中基于一維分布參數法建立了螺旋管蒸發器的動態模型, 提出合適的算法設計, 解決了計算效率的問題。在此基礎上, 針對不同輸入參數變化展開了仿真研究, 比較驗證了其相對移動邊界法的優勢。

1 模型描述

1.1 基本假設

螺線管蒸發器安裝于Li/SF6熱管反應器的外殼之中, 現將其單獨取出, 沿軸線一維展開并劃分為個控制體, 如圖1所示。

圖1 控制體劃分示意圖

在建模過程中, 對每個控制體做如下假設:

1) 由于流動而產生的沿程壓降損失是流體在螺旋管蒸發器內部的一個特性, 但在動態過程中不是關注的重點, 為減少方程數量可忽略流體在管內的壓降損失, 即不考慮動量方程;

2)工質在管內為一維流動, 忽略重力場影響;

3) 忽略工質和管壁的軸向導熱;

4) 相變流體兩相之間處于熱力學平衡狀態。

1.2 質量和能量守恒方程

對于單個控制體, 管內相變流體存在質量守恒關系

存在能量守恒關系

式中:,h,h,Q分別表示控制體內工質內能、入口比焓、出口比焓和管壁導入控制體內熱量。

1.3 壁面條件

對于單個控制體, 管壁存在能量守恒關系

式中:c,M,T,Q分別為管壁比熱容、管壁質量、管壁溫度和Li/SF6熱管反應器導入管壁的熱量。

由于螺旋管蒸發器安裝于Li/SF6反應器外殼中, 其外壁面吸收Li蒸汽凝結的潛熱, 所以外壁面邊界條件可設為勻熱流密度, 熱流密度的大小可通過Li/SF6反應器的釋熱功率核算。故外壁面換熱量

式中,q為單位長度的熱流密度。

內壁面換熱量

由于管內為湍流強制對流換熱, 且管壁與管內流體溫差較大, 故在單相區(液相區、氣相區)采用Gnielinski經驗式[12]

式中,c為溫差修正系數, 對液體

對氣體

為管內湍流流動的Darcy阻力系數, 按Filonenko經驗式

式(6)經實驗室驗證使用范圍為:=2300~ 106,=0.6~105。

對于兩相流換熱系數, 可采用V?erme Atlas經驗式[13]

式中: 下標,,,分別為飽和液體與飽和氣體;表示兩相區中氣相質量流量占總質量流量的比值。

2 模型解算

2.1 方程離散

對于微分形式的控制方程需要離散求解, 為保證求解穩定性, 采用隱格式差分。對能量守恒方程(2)進行離散, 整理成出口焓的形式

式中: 上標,1表示時刻; 控制體進出口參數的下標,1表示第個控制體的進口及出口; 控制體參數的下標表示第個控制體。

對質量守恒方程(1)進行離散, 整理成出口流量的形式

對于壁面條件, 采用顯格式離散即可。管壁能量守恒方程(3)離散得

2.2 物性方程

為使方程組閉合可求解, 還需補齊工質物性方程?，F階段的普遍做法是對工質建立物性函數庫, 在單相區直接調用, 而在兩相區則通過建立空隙率模型[14]計算得到。這使得在不同相區需使用不同方法, 導致方程數量增加, 計算效率下降。文中通過以上方式推導公式, 使兩相區計算時對物性的要求減少為只對密度和飽和溫度的求取, 再通過對美國國家標準技術研究所(national institute of standards and technology, NIST)研制開發的物性函數庫直接調用, 即可使單相區和兩相區有相同的處理方式如下

則控制體內工質質量

式中: 控制體的平均密度用式(14)求得的進出口密度的平均值即可。

2.3 解算流程

以差分格式離散微分方程解算的過程中, 不合適的算法設計易導致在空間步和時間步的計算推進過程中誤差累計, 最終導致發散, 所以設計一套合理的算法也尤為重要。文中解算流程圖如圖2所示, 算法主要分為2步, 第1步為控制體參數的預估計算, 由于采用隱格式差分, 在首次計算式會缺省部分參數, 故在首次計算時將式(11)改寫為

預估計算完控制體各參數后, 再進入第2步的循環迭代過程。在循環迭代過程中, 不停比較與上一次計算所得值的誤差, 若誤差足夠小才跳出循環, 以得到精確的控制體參數。圖2中,指上一次計算出的控制體出口質量流量, 在首次比較時,用上一時刻的控制體出口質量流量賦值;為一極小值, 大小按經驗取值, 其值大小直接影響到算法穩定性、計算精度及效率。另一影響算法穩定性的是空間步長與時間步長的取值關系, 應滿足

式中,max指工質管內流動的最大流速, 此處即為工質在螺旋管蒸發器出口處的流速。

圖2 解算流程圖

3 仿真結果與分析

基于以上數學模型, 對螺旋管蒸發器輸入參數改變時的動態過程進行仿真計算。其輸入參數分別為進入螺旋管蒸發器的工質質量流量、入口溫度和外壁面的熱流密度, 其中入口質量流量及外壁面的熱流密度分別由水泵和SF6進入反應器的流量主動調節, 入口溫則與回熱器性能相關。仿真計算時所用的穩態值及螺旋管蒸發器結構參數見表1。

表1 螺旋管蒸發器結構參數及仿真初始穩態值

3.1 入口溫階躍變化時的動態過程

3.1.1 2種模型的仿真結果

以入口溫度在0.5 s時階躍為初始溫度的95%為輸入條件, 分別得到一維分布參數模型和移動邊界模型的出口質量流量和出口比焓跟隨響應規律, 以及各相區長度的變化規律, 仿真結果見圖3~圖6。

由圖3~圖5可知, 對于一維分布參數模型入口溫度發生階躍變化后, 螺旋管蒸發器各輸出參數響應存在一個約為2.6 s的延遲, 延遲之后出口質量流量先逐漸減小后又逐漸增大, 最終穩定在與入口質量流量一致。延遲之后出口比焓先逐漸增大后又逐漸減小, 最終穩態值小于初始穩態值。延遲之后液相區長度逐漸增大, 最終穩態長度大于初始穩態值, 兩相區長度先逐漸減小后又逐漸增大, 最終穩態長度大于初始穩態值, 氣相區長度先逐漸增加后又逐漸減小, 最終穩態長度小于初始穩態值。整個動態過程持續5.5 s。

圖3 入口溫階躍時的出口質量流量

圖4 入口溫階躍時的出口比焓

圖5 一維分布參數模型各相區長度

圖6 移動邊界模型各相區長度

由圖3、圖4、圖6可知, 對于移動邊界模型入口溫度發生階躍變化后, 螺旋管蒸發器各參數發生明顯的跳躍, 跳躍后出口質量流量大于初始穩態值, 隨后經歷一個逐漸減小又逐漸增大的過程, 最終穩定在與入口質量流量一致。跳躍后出口比焓小于初始穩態值, 隨后經歷一個先逐漸增大后又逐漸減小的過程, 最終穩態值小于初始穩態值。跳躍后液相區長度小于初始穩態值, 后逐漸增大, 最終穩態長度大于初始穩態值。跳躍后兩相區長度大于初始穩態值, 隨后經歷一個先逐漸減小后又逐漸增大的過程, 最終穩態長度大于初始穩態值。跳躍后氣相區長度小于初始穩態值, 隨后經歷一個先逐漸增大后又逐漸減小的過程, 最終穩態長度小于初始穩態值。整個動態過程持續11 s。

3.1.2 仿真結果分析

對比圖3~圖6的輸入條件為入口溫度發生階躍變化后, 得到的兩模型輸出參數跟隨響應規律及各相區長度變化規律, 對兩模型結論分析對比如下。

1) 對于一維分布參數模型, 在仿真開始入口溫度發生階躍變化后, 螺旋管蒸發器各輸出參數未發生變化, 存在一個約為2.6 s的延遲, 延遲之后各輸出參數才跟隨響應。對于延遲現象及各參數發生如上響應變化的主要原因有(分析時忽略管壁的熱容效應, 假設壁面條件為均勻熱流密度): a. 由于在基礎假設中不考慮工質的軸向導熱, 故剛發生入口溫度階躍變化后, 液相區工質存在一個溫度不連續面。在溫度不連續面后的原有工質在流動到原有液體飽和位置處依舊會達到飽和, 也就是說此時液相區長度并不發生改變, 同時溫度不連續面后的原有工質參數較之前依舊保持不變, 直至溫度不連續面到達原有的液體飽和位置。這一持續時間應為工質從螺旋管蒸發器入口流至原有液體飽和位置的時間, 通過計算此時間約為2.4 s; b. 此后, 新工質到達原有液體飽和位置, 但并不足以達到飽和狀態, 故將以液體狀態繼續向前流動, 此時液相區長度逐漸增加。由于液相區密度遠大于其他2個相區, 所以螺旋管蒸發器內總工質的質量增加, 又因入口質量流量并未改變, 故此時出口質量流量減小; c. 由于出口質量流量減小, 此時還停留在螺旋管蒸發器內的原有工質流速降慢, 故只需相對初始狀態較短的兩相區長度即可加熱為飽和氣體, 且總停留時間將大于初始狀態, 故使其到達出口時具有高于之前初始狀態的溫度; d. 當原有工質完全流出螺旋管蒸發器, 出口各參數將不再變化, 動態過程結束。出口質量流量將恢復至與入口狀態相同, 出口比焓則將低于初始出口狀態。整個動態過程的時間應為工質從螺旋管蒸發器入口流至出口的時間, 通過計算約為4.2 s; e. 若考慮實際壁面條件, 不忽略管壁熱容效應, 則整個動態過程趨勢依舊如此, 但時間尺度將有所延長, 故由仿真得到的延遲時間和動態響應時間尺度上都略大于上述計算值。

綜上所述, 使用一維分布參數模型對螺旋管蒸發器發生入口溫度階躍變化進行動態仿真結果合理精確, 且該模型在仿真時未出現任何震蕩和跳躍現象, 具有很好的強壯性。

2) 使用移動邊界模型對螺旋管蒸發器入口溫度發生階躍變化進行動態仿真, 響應初期存在劇烈跳動, 并未發生延遲過程, 跳動之后的動態趨勢與一維分布參數模型所得結果基本相同, 但整體動態過程的持續時間遠長于一維分布參數模型。出現以上現象的原因主要有二, 一為入口溫度剛發生階躍變化時, 液相區中工質的溫度分布并不連續, 此時若依舊使用進出口的均值來作為液相區的平均參數將存在誤差, 這源自該模型建模時基于集中參數的思想; 其二為前文提到的仿真結果初期存在不可解釋的反沖現象。

3.2 其余輸入參數階躍變化時的動態過程

上文以入口溫度發生階躍變化, 得到的一維分布模型動態響應結果較好, 現以螺旋管蒸發器入口的工質質量流量和外壁面的熱流密度發生階躍變化為輸入條件, 得到各參數跟隨響應的動態特性, 亦是Li/SF6熱管反應器主動調節的2種方式。

3.2.1 入口質量流量階躍變化時的動態過程

以入口質量流量在0.5 s時階躍為初始質量流量的95%為輸入條件, 得到出口質量流量和出口比焓的跟隨響應規律, 仿真結果見圖7和圖8。

由圖7和圖8可知, 螺旋管蒸發器入口質量流量發生階躍變化后, 出口質量流量和出口比焓在2.5 s內并未發生較大變化, 存在一個延遲現象。2.5 s后出口質量流量逐漸減小直至與入口狀態一致, 出口比焓逐漸增大最終趨于穩定。整個動態過程持續5.5 s, 以上結果亦可反應螺旋管蒸發器入口質量流量增加時的動態特性。

圖7 入口流量階躍時的出口質量流量

圖8 入口流量階躍時的出口比焓

3.2.2 外壁面熱流密度階躍變化時的動態過程

遙感課程的考核評價由平時成績、作業成績和期末考試成績三部分組成，占比分別為20%、20%和60%。其中，平時成績主要通過考勤來考核，作業成績靠提交課后習題來考核，期末考試成績則是試卷分數。目前，對學生最終成績起決定性作用的還是期末考試成績，平時成績的認定仍停留在“上課了(學生聽與不聽，教師考核不了)，課后習題寫了(學生的創新技能考核不了)，平時成績就可以拿高分”的階段。這種考核評價方法沒有考慮到學生的學習過程，無法真實地體現學生對知識的消化程度和對實踐技能的掌握程度，不能反映學生的真實水平，也不能培養學生解決實際問題的能力。

以外壁面熱流密度在0.5 s時階躍為初始熱流密度的95%為輸入條件, 得到出口質量流量和出口比焓的跟隨響應規律, 仿真結果見圖9和圖10。

圖9 熱流密度階躍時的出口質量流量

圖10 熱流密度階躍時的出口比焓

由圖9和圖10可知, 螺旋管蒸發器外壁面熱流密度發生階躍變化后, 出口質量流量立即突變減小, 在接下來的2.5 s內, 出口質量流量和出口比焓并未發生較大變化, 存在一個延遲現象。2.5 s后, 出口質量流量逐漸增大直至與入口狀態一致, 出口比焓逐漸減小最終趨于穩定。整個動態過程持續5.5 s, 以上結果亦可反應螺旋管蒸發器外壁面熱流密度增加時的動態特性。

3.2.3 兩模型仿真時間對比

使用同一計算機和相同的時間步長進行15 s的仿真過程中, 一維分布參數模型共用時59 s計算出結果, 而移動邊界模型用時48 s?？梢娢闹兴龅囊痪S分布參數模型動態仿真過程具有與移動邊界模型接近的運算效率, 可用于系統級的動態仿真之中。

4 結論

文中基于一維分布參數法針對螺旋管蒸發器建立了相應模型, 對輸入參數為入口質量流量、入口溫度和外壁面熱流密度發生階躍變化進行了仿真計算, 并對比移動邊界模型的結果, 得到如下結論。

1) 建立的一維分布參數模型能明確描述螺旋管蒸發器動態初期的延遲現象和整體的動態趨勢及時間尺度, 能較為準確地反映螺旋管蒸發器在輸入參數改變時的動態響應過程。并且在輸入參數劇烈變化時模型計算依舊穩定, 未出現任何跳躍、震蕩現象, 具有很好的強壯性。

2) 運用移動邊界模型對螺旋管蒸發器輸入參數改變時進行動態仿真初期存在劇烈的跳躍, 跳動之后的動態響應趨勢與一維分布參數模型相同, 但時間尺度遠大于一維分布參數模型。其原因主要源自建模時基于集中參數的思想和仿真過程中出現的不可解釋的反沖現象。

3) 通過合適的建模方式和解算流程, 一維分布參數模型能得到較高的計算效率, 在同等情況下僅落后移動邊界模型18.6%, 使其能用于實際的系統級動態仿真, 并為后期控制系統設計打下基礎。

[1] 查志武. 魚雷熱動力技術[M]. 北京: 國防工業出版社, 2006.

[2] Blakeslee T R. Laminar Free-convective Combustion of a Liquid Metal from a Wick[D]. U.S: Pennsylvania State University, 1977.

[3] Zhang H Q, Zhou L X, Chan C K. A Double--, and Multi-δ, PDFModel for Turbulent Gas-liquid Reacting Flows [J]. Fuel, 2002, 81(6): 817-827.

[4] 李侃侃, 程惠爾, 臧家亮. 用均質流模型研究SF_6/Li燃燒體系[J]. 上海交通大學學報, 1999, 33(8): 1017- 1019.Li Kan-kan, Cheng Hui-er, Zang Jia-liang. Study on SF6/Li Combustion System Using the Model of Homo- geneous Flow[J]. Journal of Shanghai Jiaotong University, 1999, 33(8): 1017-1019.

[5] 鄭邯勇, 卜建杰. 六氟化硫在熔融鋰中的浸沒噴射反應過程[J]. 化工學報, 1996, 47(6): 656-662.Zheng Han-yong, Bu Jian-jie. The Submerged Jet Reac- tion Process of Sulfur Hexafluoride into Molten Lithi- um[J]. CIESC Journal, 1996, 47(6): 656-662.

[6] He X, Liu S, Asada H. Modeling of Vapor Compression Cycles for Advanced Controls in HVAC systems[C]//Am- erican Control Conference. U.S: IEEE, 1995.

[7] Zhang W J, Zhang C L. 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.

[8] 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.

[9] 白杰. 無人水下航行器新型熱電聯合閉式循環熱動力系統研究[D]. 西安: 西北工業大學, 2016.

[10] Jia X, Tso C P, Chia P K, et al. A Distributed Model for Prediction of the Transient Response of an Evaporator[J]. International Journal of Refrigeration, 1995, 18(5): 336- 342.

[11] Wang C C. A Numerical Method for Thermally Non- equilibrium Condensing Flow in a Double-pipe Conden- ser[J]. Applied Thermal Engineering, 1997, 17(7): 647- 660.

[12] 楊世銘, 陶文銓. 傳熱學[M]. 第4版. 北京: 高等教育出版社, 2006.

[13] Jensen J M. Dynamic Modeling of Thermo-Fluid Systems with Focus on Evaporators for Refrigeration[D]. Denm- ark: Technical University of Denmark, 2003.

[14] Rice C K. The Effect of Void Fraction Correlation and Heat Flux Assumption on Refrigerant Charge Inventory Pr- edictions[J]. Ashrae Transactions, 1987, 93(1): 341-367.

(責任編輯: 楊力軍)

Dynamic Modeling of Spiral Tube Evaporator Based on One-Dimensional Distributed Parameter Method

MIN Zhi-yong, DANG Jian-jun, SU Hao, ZHANG Jia-nan, LI Dai-jin

(School of Marine Science and Technology, Northwestern Polytechnical University, Xi′an 710072, China)

The metal fuel with higher energy density has become an important role in development of underwater thermal power systems. In this paper, the one-dimensional distributed parameter method is employed to establish a dynamic model of the spiral tube evaporator in a Li/SF6heat pipe reactor. Simulations are performed by changing different input parameters, and the results are compared with that of the moving boundary method. It is indicated that the one-dimen- sional distributed parameter model can better represent the real system, and is able to reflect the dynamic response of the system more accurately. The computational efficiency approaches that of the moving boundary model if suitable formula derivation and proper algorithm are adopted. This model is robust with respect to sudden change of the input load. This research may provide the basis for design, control and real-time simulation of underwater thermal power systems.

underwater thermal power system; spiral tube evaporator; one-dimensional distributed parameter method; moving boundary method; dynamic modeling

閔智勇,黨建軍,蘇浩,等.基于一維分布參數法的螺旋管蒸發器動態建模[J].水下無人系統學報,2017,25(4):344-350.

TJ630.1; TL353.1

A

2096-3920(2017)04-0344-07

10.11993/j.issn.2096-3920.2017.04.007

2017-06-21;

2017-07-14.

國家自然科學基金面上項目(51679202), 西北工業大學研究生創意創新種子基金資助(Z2017079).

閔智勇(1993-), 男, 在讀碩士, 主要研究方向為水下熱動力系統仿真及控制.