王海濤,羅靖,張麗敏,孫巖雷,陳遠(yuǎn)航,常德鵬,陳燕燕,余國瑤,羅二倉

(1.中國科學(xué)院理化技術(shù)研究所低溫工程學(xué)重點(diǎn)實(shí)驗(yàn)室,100190,北京;2.中國科學(xué)院大學(xué),100049,北京)

熱聲熱機(jī)也稱回?zé)崾綗釞C(jī),因具有高效率、高可靠性、長壽命等優(yōu)勢(shì)而逐漸受到重視[1-3]。換熱器是實(shí)現(xiàn)發(fā)動(dòng)機(jī)內(nèi)部工作氣體與外界熱量交換的核心部件,回?zé)崞魑挥诟?、低溫?fù)Q熱器之間,實(shí)現(xiàn)熱能與聲能的相互轉(zhuǎn)換。高、低溫?fù)Q熱器與回?zé)崞鞯炔考o密相連,工質(zhì)在相鄰部件間往復(fù)穿梭。如圖1所示,在熱聲發(fā)動(dòng)機(jī)的換熱器中,絕大部分工質(zhì)在各交界面的突變截面處往復(fù)流動(dòng),顯著的時(shí)變慣性與“進(jìn)、出口”效應(yīng)使得熱聲熱機(jī)換熱器的特性必然與定常流動(dòng)換熱器存在顯著不同。

圖1 熱聲發(fā)動(dòng)機(jī)核心部件示意圖

由于物理過程本身的復(fù)雜性與多樣性,至今對(duì)交變流動(dòng)及其換熱規(guī)律的認(rèn)識(shí)與定常流相比仍然相差甚遠(yuǎn)。理論分析工作主要基于線性或弱非線性熱聲理論,針對(duì)層流充分發(fā)展換熱器部件進(jìn)行分析[4-7]。實(shí)驗(yàn)與數(shù)值分析研究方面,代表性工作主要集中在20世紀(jì)80—90年代對(duì)斯特林熱機(jī)的加熱器或冷卻器的研究[8-12],以及近些年來的二維簡單板疊式換熱器的仿真與可視化[13-17]、三維瞬態(tài)仿真研究[18-20]等。前一階段的研究主要是借鑒定常流的分析方法,嘗試采用阻力系數(shù)與換熱系數(shù)來描述換熱器的特征,但所得到的表達(dá)式往往相差較大,缺乏普適性,這其中的緣由也鮮少有分析。后一階段則主要針對(duì)簡單換熱器進(jìn)行CFD二維與三維模擬以及流場(chǎng)可視化測(cè)量,考察不同聲場(chǎng)強(qiáng)度下?lián)Q熱器內(nèi)的瞬態(tài)流動(dòng)、換熱分布特征,以及換熱器端部的渦流問題等。這種細(xì)致分析的結(jié)果顯示出交變流動(dòng)換熱器與定常流換熱器有著顯著的不同,首先是熱流密度的分布與聲場(chǎng)分布及系統(tǒng)其他部件之間的耦合特性非常強(qiáng);另外,顯示出換熱器的流動(dòng)與換熱特征主要受進(jìn)出口效應(yīng)影響。但是,這類研究未能在流場(chǎng)細(xì)節(jié)之外總結(jié)出規(guī)律,且受限于研究的具體熱機(jī)形式及實(shí)驗(yàn)條件(針對(duì)的是簡單換熱器結(jié)構(gòu)、單一聲場(chǎng)條件以及低功率密度下的研究),分析結(jié)果的系統(tǒng)性、普適性和實(shí)用性都顯得不足。

隨著大功率熱聲熱機(jī)發(fā)展的需求日益迫切,實(shí)際熱源耦合的換熱器設(shè)計(jì)越來越成為關(guān)鍵,然而對(duì)熱機(jī)內(nèi)部交變流動(dòng)換熱器特性的理解尚難以滿足發(fā)展需求。因此,有必要對(duì)熱聲熱機(jī)換熱器內(nèi)部的換熱和阻力特征展開細(xì)致研究,明確換熱器內(nèi)部換熱與流動(dòng)損失特性,系統(tǒng)分析并總結(jié)特性規(guī)律,有效指導(dǎo)熱聲熱機(jī)的工程設(shè)計(jì)與應(yīng)用。本文以一臺(tái)新型10 kW熱管耦合自由活塞斯特林發(fā)電機(jī)的發(fā)動(dòng)機(jī)核心單元為分析對(duì)象,分別采用回?zé)崾綗釞C(jī)通用設(shè)計(jì)準(zhǔn)一維軟件Sage[21-22]和商用CFD計(jì)算軟件FLUENT[23],對(duì)高、低溫?fù)Q熱器內(nèi)的流動(dòng)與換熱特征進(jìn)行詳細(xì)分析?；谕?fù)涞刃Ｐ?對(duì)比分析兩種計(jì)算結(jié)果,提出適用于Sage計(jì)算的換熱器修正模型,提升Sage模型的計(jì)算準(zhǔn)確性,指導(dǎo)熱聲熱機(jī)的工程設(shè)計(jì)。

1 熱管耦合發(fā)電機(jī)整體設(shè)計(jì)

1.1 熱管耦合自由活塞斯特林發(fā)電機(jī)結(jié)構(gòu)

熱管耦合自由活塞斯特林發(fā)電機(jī)由實(shí)現(xiàn)熱聲轉(zhuǎn)換的發(fā)動(dòng)機(jī)與聲電轉(zhuǎn)換的直線電機(jī)兩部分組成。發(fā)動(dòng)機(jī)側(cè)需要實(shí)現(xiàn)熱管與發(fā)動(dòng)機(jī)加熱器的結(jié)構(gòu)與熱耦合,熱管有柱形直熱管(亦可大角度彎折)與異形熱管兩種,如圖2所示。柱形直熱管更為成熟,與自由活塞斯特林的耦合方式包括熱頭橫截面插入與軸向插入式兩種。常規(guī)耦合受限于發(fā)動(dòng)機(jī)熱頭在橫截面與軸向長度上的過度緊湊,使得熱管冷凝段長度與發(fā)動(dòng)機(jī)的尺寸適配性較差,只能用于低功率匹配,對(duì)于高功率、緊湊型發(fā)動(dòng)機(jī)結(jié)構(gòu)的應(yīng)用存在限制。異形集成式熱管加熱器能實(shí)現(xiàn)發(fā)動(dòng)機(jī)與熱管之間的高效耦合,發(fā)電機(jī)功率可達(dá)幾十千瓦。熱管結(jié)構(gòu)與發(fā)動(dòng)機(jī)的承壓壁集成,異形熱管具有顯著結(jié)構(gòu)依賴特點(diǎn),通用性差,結(jié)構(gòu)復(fù)雜,工藝要求高,目前尚不成熟。

針對(duì)以上問題,本文設(shè)計(jì)出一種高效軸向直熱管自由活塞斯特林發(fā)電機(jī),通過降低加熱器孔隙率與降低系統(tǒng)工作頻率,提高加熱器內(nèi)工質(zhì)位移幅度,從而有效增長加熱器長度,即增長熱管冷凝段長度。本文所研究的發(fā)電機(jī)結(jié)構(gòu)優(yōu)化設(shè)計(jì)基于Sage軟件完成,如圖3所示,為基于熱力-電磁設(shè)計(jì)完成的整機(jī)基本結(jié)構(gòu)示意圖。直熱管插入如圖3(b)所示的套筒內(nèi),氣體在銅塊的矩形窄縫中進(jìn)行往復(fù)流動(dòng)換熱,同時(shí)為實(shí)現(xiàn)熱管加熱器與發(fā)電機(jī)的高效耦合,負(fù)載發(fā)電機(jī)部分如圖3(a)所示,采用對(duì)置電機(jī)的形式,實(shí)現(xiàn)單側(cè)發(fā)電機(jī)熱耦合的同時(shí)減小系統(tǒng)振動(dòng)。

(a)發(fā)動(dòng)機(jī)整體結(jié)構(gòu) (b)熱管加熱器局部

經(jīng)設(shè)計(jì)優(yōu)化后的換熱器結(jié)構(gòu)參數(shù):翅片式高溫加熱器長度為250 mm,流道寬度為1.3 mm,流道高度為10 mm,翅片平均厚度為7.986 mm,流道數(shù)為108,孔隙率為0.046;回?zé)崞鏖L度為52 mm,采用絲綿填充,孔隙率為0.897;室溫?fù)Q熱器為管殼式,長度為70 mm,管徑為2 mm,孔隙率為0.08,流道數(shù)為768。

1.2 整機(jī)性能參數(shù)

熱管耦合發(fā)電機(jī)整機(jī)優(yōu)化后的額定工況性能參數(shù)如表1所示。額定設(shè)計(jì)發(fā)電功率為10.77 kW,對(duì)應(yīng)加熱器吸熱量為30 kW,冷卻器排熱量約為19 kW。

表1 熱管耦合自由活塞斯特林發(fā)電機(jī)性能參數(shù)

發(fā)動(dòng)機(jī)側(cè)冷卻器、回?zé)崞髋c加熱器核心段的壓力、體積流率波動(dòng)幅值分布如圖4(a)所示,聲功、聲阻抗相位分布如圖4(b)所示。為適應(yīng)熱管冷凝段長度,加熱器長度為250 mm,可見加熱器內(nèi)明顯的壓力幅值下降與聲功降低,即加長的換熱器明顯增加了聲功損失。

(a)發(fā)動(dòng)機(jī)內(nèi)壓力波動(dòng)與體積流率分布

基于以上熱管耦合發(fā)電機(jī)結(jié)構(gòu)與發(fā)動(dòng)機(jī)核心部件額定工況,本文以冷卻器、回?zé)崞髋c加熱器組成的核心單元為研究對(duì)象,耦合熱聲轉(zhuǎn)換,重點(diǎn)分析加熱器與冷卻器內(nèi)的流動(dòng)與換熱特性。

2 換熱器計(jì)算模型

2.1 換熱器拓?fù)涞刃?/h3>

為實(shí)現(xiàn)對(duì)換熱器內(nèi)流動(dòng)與換熱特征的系統(tǒng)分析,同時(shí)簡化計(jì)算過程,將三維流體計(jì)算區(qū)域通過拓?fù)涞刃Ш喕癁槎S平面結(jié)構(gòu)。加熱器與冷卻器均為窄縫結(jié)構(gòu),相鄰的部件為回?zé)崞髋c空腔(加熱器側(cè)連接膨脹腔,冷卻器側(cè)連接壓縮腔),回?zé)崞鳛槎嗫捉橘|(zhì)結(jié)構(gòu),具有顯著的均勻化流動(dòng)效果。因此,回?zé)崞髋c冷卻器的拓?fù)鋵?duì)應(yīng)性即可解耦,即將加熱器與冷卻器窄縫簡化為二維結(jié)構(gòu)時(shí),與回?zé)崞鞯慕Y(jié)構(gòu)匹配只需要滿足界面上的孔隙率對(duì)應(yīng)即可。由于換熱器孔隙率相對(duì)于回?zé)崞髋c空腔足夠小,因此忽略換熱器窄縫三維分布對(duì)回?zé)崞髋c空腔內(nèi)射流與二維射流的差異。以狹長窄縫二維等效為基礎(chǔ),以孔隙率最小的加熱器一條完整窄縫為最小計(jì)算單元,以窄縫寬度為基準(zhǔn)。冷卻器原設(shè)計(jì)為管束,二維近似通過水力直徑、換熱面積及孔隙率相近原則進(jìn)行等效,等效后的窄縫數(shù)量取相對(duì)加熱器鄰近整數(shù),適當(dāng)調(diào)整水力直徑以匹配相對(duì)孔隙率,回?zé)崞鲗挾扰c空腔寬度以相對(duì)換熱器孔隙率確定,加熱器、回?zé)崞?、冷卻器的長度保持不變,空腔長度則以相對(duì)容積確定。從而將三維結(jié)構(gòu)等效為二維平面結(jié)構(gòu)。壓縮腔與膨脹腔容積根據(jù)CFD動(dòng)活塞運(yùn)動(dòng)要求進(jìn)行容積擴(kuò)展,以壓縮腔與冷卻器界面、膨脹腔與加熱器界面上的壓力、體積流率與整機(jī)完全相同為原則,將計(jì)算邊界外延。

圖5為按照拓?fù)涞刃г瓌t獲得的自由活塞斯特林發(fā)電機(jī)二維模型,結(jié)構(gòu)參數(shù)如表2所示?；谠摱S模型,可同時(shí)采用sage與CFD進(jìn)行仿真計(jì)算,sage一維模型計(jì)算結(jié)果即可以與原整機(jī)設(shè)計(jì)計(jì)算結(jié)果進(jìn)行對(duì)比校驗(yàn),以檢驗(yàn)拓?fù)浜喕^程的合理性,同時(shí)也可為CFD局部仿真提供邊界條件?；贑FD仿真對(duì)換熱器時(shí)變換熱特征進(jìn)行詳細(xì)分析,并與sage一維模型計(jì)算結(jié)果進(jìn)行對(duì)比,考察sage一維模型計(jì)算的有效性。

表2 等效換熱核心單元結(jié)構(gòu)參數(shù)

圖5 發(fā)動(dòng)機(jī)核心單元拓?fù)涞刃Ｐ?/p>

圖5所示計(jì)算模型中,高溫加熱器的壁面溫度為823 K,低溫冷卻器的壁面溫度為333 K,其它壁面邊界設(shè)置均為絕熱。軸向計(jì)算邊界為避免開口邊界帶來的入流焓流不確定性,采用雙活塞邊界進(jìn)行計(jì)算。活塞邊界即設(shè)定往復(fù)振蕩虛擬活塞邊界面,模擬活塞對(duì)熱聲核心單元的往復(fù)壓縮、膨脹效應(yīng),能量從邊界上以時(shí)均聲功形式進(jìn)出系統(tǒng),活塞邊界可設(shè)置絕熱,無不確定性能量,因此對(duì)于系統(tǒng)內(nèi)部的能量分布分析精度更高。Sage一維模型可采用虛擬活塞邊界進(jìn)行計(jì)算;CFD仿真中,兩端虛擬活塞邊界則必須采用動(dòng)網(wǎng)格,為保證活塞相對(duì)平衡位置不變,還需要進(jìn)行模擬充氣與活塞同步的初始過程模擬?；钊胶?即可加載壓縮活塞與膨脹活塞運(yùn)動(dòng)的速度控制,活塞速度隨時(shí)間的變化函數(shù)如下

upiston(t)=ωxpiston1cos(ωt+θpiston)

(1)

式中:xpiston1、θpiston分別表示活塞面的位移幅值與相位;ω是角頻率。

2.2 Sage一維模型

Sage是由美國Gedeon等于1995年根據(jù)MS-DOS編寫開發(fā)的一款商業(yè)計(jì)算軟件,主要是針對(duì)回?zé)崾綗釞C(jī)的通用設(shè)計(jì)計(jì)算軟件。其采用圖形化界面將與實(shí)際物理組件相對(duì)應(yīng)的氣體流動(dòng)、傳熱和其他建模細(xì)節(jié)封裝在一些特定的模型組件中,用戶通過將封裝好的組件按照特定的需求進(jìn)行相互連接,即可對(duì)復(fù)雜的回?zé)崾綗釞C(jī)系統(tǒng)建立仿真模型,利用質(zhì)量、動(dòng)量和能量守恒方程、本構(gòu)方程將各部件的參數(shù)進(jìn)行有機(jī)耦合,最終實(shí)現(xiàn)整機(jī)的求解與多參數(shù)協(xié)同優(yōu)化。圖6給出圖5所對(duì)應(yīng)的Sage一維計(jì)算模型。

ρstdy—系統(tǒng)載荷;Pphsr—體積位移;mGt—進(jìn)出口的質(zhì)量流;Qstdy—線性熱源;T0—低溫?zé)嵩?Th—高溫?zé)嵩础?/p>

回?zé)崞髂Ｐ瓦x擇軟件內(nèi)嵌的絲綿多孔介質(zhì)模型,其阻力的本構(gòu)方程如下

f=a1/Re+a2Rea3

(2)

式中:a1=25.7α+79.8;a2=0.146α+3.76;a3=-0.002 83α-0.074 8;α=ε/(1-ε);ε為孔隙率。

換熱特性的本構(gòu)方程如下

Nu=1+b1Peb2

(3)

式中:b1=0.186α;b2=0.55;Pe為貝克來數(shù),即雷諾數(shù)Re與普朗特?cái)?shù)Pr的乘積。

等效模型中的換熱器均為平板模型,采用軟件內(nèi)嵌的板式換熱器模型,基于線性熱聲理論平板模型獲得板式換熱器模型的阻力與換熱特征層流,湍流則借用湍流平板穩(wěn)態(tài)流動(dòng)特性,本構(gòu)方程如下

f=0.11(ε/dh+68/Re)0.25

(4)

Nu=0.035Re0.75Pr0.33

(5)

式中:dh為絲綿等效水力直徑。

Sage一維模型中可通過設(shè)置額外的阻力系數(shù),對(duì)模塊進(jìn)行整體阻力修正,阻力修正通過定義以下壓力梯度表達(dá)式中的K實(shí)現(xiàn)。即計(jì)算模塊的壓降由模塊特征阻力系數(shù)f與附加可用戶自定義的阻力系數(shù)K決定。

在流動(dòng)方向上的壓降表達(dá)式如下

(6)

式中:L為模型長度;ρu|u|/2為流體動(dòng)能。

2.3 CFD模型

CFD仿真使用ANSYS FLUENT 2021 R1版本。采用活塞邊界動(dòng)網(wǎng)格進(jìn)行瞬態(tài)計(jì)算,發(fā)動(dòng)機(jī)內(nèi)氦氣的交變流動(dòng)計(jì)算屬于小馬赫數(shù)可壓縮流動(dòng)。將氦氣設(shè)置為理想氣體,采用基于壓力法的求解器計(jì)算,選用PISO 算法,壓力使用PRESTO!離散格式,其他均采用二階迎風(fēng)離散格式,時(shí)間離散采用二階時(shí)間差分。換熱器內(nèi)流速較高,以幅值雷諾數(shù)作為判斷基準(zhǔn),加熱器雷諾數(shù)幅值變化范圍在3 000以上,明顯高于2 300,因此可采用k-ε湍流模型。模型中計(jì)入氦氣物性隨溫度的變化,物性公式列于表3?；?zé)崞鞑捎脽崞胶舛嗫捉橘|(zhì)模型,只考慮回?zé)崞鲀?nèi)流阻模型。阻力系數(shù)根據(jù)Sage模型中所用的阻力系數(shù)公式,通過局部擬合轉(zhuǎn)換為黏性阻力系數(shù)與慣性阻力系數(shù)表征的阻力特性,獲得FLUENT中多孔介質(zhì)模型設(shè)置參數(shù)1/K(黏性阻力系數(shù))與C2(慣性阻力系數(shù))。

表3 氦氣物性參數(shù)

2.4 CFD二維模型網(wǎng)格與時(shí)間尺度無關(guān)性檢驗(yàn)

基于額定工況,進(jìn)行加密網(wǎng)格與時(shí)間尺度的計(jì)算,同時(shí)對(duì)壓力、體積流率、溫度、聲功、換熱量的波動(dòng)與時(shí)均值進(jìn)行監(jiān)測(cè),波動(dòng)量幅值、相位及時(shí)均量均通過FFT分析獲得。圖7與圖8給出了典型參數(shù)的影響。其中網(wǎng)格尺度影響相對(duì)較小,網(wǎng)格量大于19 800后,對(duì)計(jì)算結(jié)果影響可忽略,因此后續(xù)計(jì)算網(wǎng)格均采用19 800這套網(wǎng)格模型。

圖7 壓力波動(dòng)幅度延程分布隨網(wǎng)格數(shù)的變化

圖8 不同時(shí)間步長下時(shí)均能量守恒性隨計(jì)算時(shí)間的變化

時(shí)間尺度對(duì)波動(dòng)參數(shù)幅度的影響較小,但是對(duì)時(shí)均能量的影響較大。由于計(jì)算區(qū)域?yàn)榉忾]系統(tǒng),活塞邊界的時(shí)均聲功差異應(yīng)該等于高、低溫?fù)Q熱器的時(shí)均換熱量差異,但CFD二維模型的時(shí)均能量計(jì)算發(fā)現(xiàn)存在明顯的能量不平衡,且該能量不平衡性直接與時(shí)間尺度相關(guān)。如圖8中黑色線代表的每周期計(jì)算500個(gè)點(diǎn)時(shí),系統(tǒng)有600 W左右的不平衡;每周期計(jì)算達(dá)到5 000個(gè)點(diǎn)時(shí),系統(tǒng)能量不平衡降低到170 W左右。繼續(xù)減小時(shí)間步長,會(huì)導(dǎo)致計(jì)算量過大,相較于10 kW量級(jí)的系統(tǒng)典型時(shí)均能流,計(jì)算誤差可接受,本文選擇一個(gè)周期計(jì)算5 000步。對(duì)小時(shí)間步長的要求意味著熱聲系統(tǒng)中時(shí)均能量的準(zhǔn)確求解較為困難,這主要是因?yàn)槟芰康闹饕卣鳛椴▌?dòng)量,而時(shí)均值為多參數(shù)耦合值,相較于波動(dòng)量為高階小量,且與參數(shù)間的相位關(guān)系直接相關(guān),因而計(jì)算精度對(duì)時(shí)間尺度要求非常高。

2.5 額定工況Sage模型與CFD模型計(jì)算對(duì)比

表4給出了基于圖5的Sage模型與CFD模型計(jì)算結(jié)果對(duì)比。其中:Pcomp表示壓縮腔入口壓力一階幅值;Pexp1表示膨脹腔1入口壓力一階幅值;PVchxout表示冷卻器出口截面時(shí)均聲功;PVhhxin表示加熱器入口截面時(shí)均聲功;Qchx表示冷卻器壁面時(shí)均換熱量;Qhhx表示加熱器壁面時(shí)均換熱量。兩種模型計(jì)算均采用相同活塞邊界,因此壓縮腔與膨脹腔活塞界面處的體積流率完全相同。二者計(jì)算的時(shí)均能量相當(dāng),差異在10%附近,但壓力差異相對(duì)較大,即內(nèi)部阻力特性差異較大。

1、20世紀(jì)50年代初的探索期是箏樂演奏技術(shù)發(fā)展的一個(gè)分水嶺，打破了“右手演奏，左手和聲”的觀念。這個(gè)時(shí)期的代表作《慶豐年》。

表4 兩種模型相同截面計(jì)算參數(shù)的對(duì)比

3 換熱器特性分析

為詳細(xì)分析發(fā)動(dòng)機(jī)換熱器的換熱特性,基于固定結(jié)構(gòu)參數(shù)與額定工況,保持其他參數(shù)不變,進(jìn)行不同位移幅度與不同聲場(chǎng)相位下的換熱器內(nèi)流動(dòng)與換熱特性分析。其中位移幅值的變化是通過等比例改變活塞速度來實(shí)現(xiàn)的,從額定工況的0.1倍變化到2倍,改變位移幅度時(shí),計(jì)算區(qū)域內(nèi)的聲場(chǎng)相位分布變化非常小。聲場(chǎng)相位的變化通過保持額定工況位移幅度不變、改變兩個(gè)活塞之間的相位差(在0～2π范圍內(nèi)變化)實(shí)現(xiàn)。

3.1 換熱器內(nèi)換熱分布特性

如圖9所示的額定工況下?lián)Q熱分布對(duì)比可知,CFD二維模型與Sage一維模型的計(jì)算對(duì)時(shí)均溫差與時(shí)均熱流密度qw分布的捕捉基本相似,數(shù)值相近,說明Sage一維模型可以合理地獲得換熱器內(nèi)時(shí)均換熱分布與集中特征。但是,與CFD二維模型計(jì)算結(jié)果仍然存在明顯差異,該差異與溫度場(chǎng)的局部二維特性有關(guān),也與流阻損失產(chǎn)生的場(chǎng)能量分布差異有關(guān)。計(jì)算聲場(chǎng)下,換熱器表現(xiàn)為靠近回?zé)崞饕粋?cè)換熱更為顯著,由于為等壁溫條件,因此對(duì)應(yīng)的時(shí)均溫差也更大。

(a)額定工況下?lián)Q熱器內(nèi)時(shí)均溫差分布

圖10中不同位移幅度下的時(shí)均熱流密度分布變化表明,隨著換熱器內(nèi)位移幅度低于額定工況,換熱器內(nèi)出現(xiàn)無效換熱區(qū)間(即時(shí)間換熱為0),超過額定工況后,換熱器換熱強(qiáng)度分布出現(xiàn)近等比例抬升。Sage一維模型計(jì)算與CFD二維模型計(jì)算結(jié)果結(jié)論一致,說明換熱器有效長度必須與其內(nèi)部工質(zhì)位移幅度相匹配,位移過小,換熱器有效換熱區(qū)域變小,印證了熱聲系統(tǒng)中的換熱器進(jìn)、出口效應(yīng)本質(zhì)與本征緊湊性特征。

圖10 不同位移幅度下時(shí)均熱流密度分布的CFD二維模型計(jì)算結(jié)果

圖11將時(shí)均熱流密度與時(shí)均氣固溫差代入努塞特?cái)?shù)計(jì)算式,獲得時(shí)均等效Nu(Nu0),冷卻器中的Nu0顯著高于加熱器,且延程分布特性更為顯著,對(duì)速度也更為敏感。這進(jìn)一步說明了熱聲發(fā)動(dòng)機(jī)中換熱器換熱特征具有顯著的局部特征。

圖11 不同位移幅度下Nu0的CFD二維模型計(jì)算結(jié)果

鑒于Sage一維模型能較好地捕捉換熱器內(nèi)的換熱特征,且計(jì)算速度顯著快于CFD二維模型,因此以下采用Sage一維模型計(jì)算不同聲場(chǎng)相位時(shí)的換熱器內(nèi)換熱特性。如圖12所示,調(diào)節(jié)計(jì)算邊界上的活塞相位差(Δθ),得到的各工況下聲場(chǎng)相位分布圖。圖13給出相應(yīng)的時(shí)均聲功分布情況,可見核心單元內(nèi)的能流分布發(fā)生了顯著變化。

圖12 不同邊界活塞相位差下的聲阻抗角分布

圖13 不同邊界活塞相位差下的時(shí)均聲功分布

由圖14可知,雖然加熱器與冷卻器外壁面溫度保持不變,即加熱器外壁面為823 K,冷卻器外壁面為333 K,但聲場(chǎng)的變化不僅改變了熱流密度分布特性,還決定了換熱器是吸熱還是放熱。因此,熱聲系統(tǒng)中的換熱器換熱特性是聲場(chǎng)能流驅(qū)動(dòng)的,由時(shí)均焓流決定時(shí)均換熱,而不是溫差驅(qū)動(dòng)換熱。時(shí)均溫差是時(shí)均能流變化產(chǎn)生時(shí)均換熱的結(jié)果而不是換熱的驅(qū)動(dòng)力,這與穩(wěn)態(tài)流動(dòng)完全不同。

圖14 不同邊界活塞相位差下的時(shí)均熱流密度分布

3.2 阻力特性分析

圖15給出了額定工況下冷卻器、回?zé)崞鳌⒓訜崞髦袎毫?、體積流率、聲阻抗相位角以及時(shí)均聲功的分布對(duì)比。其中最大的差異是壓力波動(dòng)在回?zé)崞髋c換熱器的交界面上的壓力差異。這是因?yàn)镾age作為一維模型無法直接獲得復(fù)雜流動(dòng)變化帶來的局部損失,且已有換熱器模塊不能針對(duì)部件間的相對(duì)結(jié)構(gòu)特征進(jìn)行自動(dòng)局部阻力調(diào)整,而是需要人為設(shè)定相關(guān)損失。CFD二維模型仿真可以獲得突變截面流場(chǎng)特性的變化,包括射流的影響。從計(jì)算結(jié)果可以看出,由于加熱器孔隙率非常小,在回?zé)崞髋c加熱器界面上產(chǎn)生了顯著的一階波動(dòng)壓力損失,時(shí)均能量也可看到明顯的聲功損失。換熱器與兩側(cè)腔體的連接處損失則相對(duì)小得多。這是因?yàn)?加熱器的高速射流在回?zé)崞鬟@一高阻力部件中產(chǎn)生的損失要比射流在空腔中顯著得多。這是換熱器計(jì)算模型中Sage一維模型與CFD二維模型的最大差異來源。

(a)壓力

(a)壓力

3.3 Sage模型阻力修正

基于CFD模型計(jì)算結(jié)果,可對(duì)換熱器與回?zé)崞鞯耐蛔兘缑嫣幍牧髯钃p失進(jìn)行特征擬合,反代入Sage模型計(jì)算模型進(jìn)行修正,將修正結(jié)果再次與CFD模型計(jì)算結(jié)果進(jìn)行對(duì)比,即可判斷修正模型的適用性。另一方面,基于穩(wěn)態(tài)突變界面局部損失系數(shù)與Swift計(jì)算方法亦可對(duì)突變界面特性進(jìn)行修正。以下將針對(duì)這兩類方法進(jìn)行對(duì)比驗(yàn)證。

由Swift方法對(duì)熱聲系統(tǒng)中噴射泵相關(guān)理論分析[24]可知,氣體往復(fù)穿梭變截面時(shí),半個(gè)周期經(jīng)歷流道擴(kuò)張,半個(gè)周期經(jīng)歷流道收縮,變截面處的局部阻力系數(shù)K應(yīng)具有時(shí)間依賴性,即為時(shí)間的函數(shù),因此

(7)

使用正弦速度

U(t)=|U1|sinωt

(8)

將式(8)代入式(7),可求得平均壓降為

(9)

則其微分形式為

(10)

根據(jù)穩(wěn)態(tài)流動(dòng)中面積比與管道截面突變局部阻力系數(shù)表[25],計(jì)算模型面積比,即可獲得交變流動(dòng)界面等效阻力系數(shù)值。

根據(jù)CFD模型計(jì)算所得壓力波動(dòng)分布進(jìn)行等效阻力系數(shù)修正嘗試了兩種方法:一種是針對(duì)波動(dòng)壓力幅值與速度幅值進(jìn)行修正;另一種是針對(duì)波動(dòng)壓力與速度的有效值進(jìn)行修正。其中波動(dòng)壓力幅值修正通過下式計(jì)算出CFD二維模型的波動(dòng)幅值壓降與速度幅值之間的關(guān)系,再計(jì)算出總阻力系數(shù)Ktotal

(11)

(12)

對(duì)于有效值計(jì)算法則,根據(jù)下式計(jì)算CFD二維模型的換熱器壓力與速度瞬時(shí)值,再計(jì)算出總阻力系數(shù)Ktotal

(13)

通過以上方法,針對(duì)冷卻器與加熱器在不同速度幅值下的總阻力系數(shù)進(jìn)行計(jì)算。Swift計(jì)算方法獲得的阻力系數(shù)與速度無關(guān),CFD模型數(shù)據(jù)處理獲得兩種阻力系數(shù)是速度的函數(shù)。將以上3種修正方法進(jìn)行對(duì)比發(fā)現(xiàn),采用有效壓降法獲得的修正最為有效?；谟行毫π拚椒ㄊ沟肧age模型波動(dòng)壓降與聲功損耗均高于CFD模型計(jì)算結(jié)果,但差異存在比例關(guān)系。因此,基于有效壓力的修正方法需再引入一個(gè)阻力系數(shù)修正因子

φ≈0.7

(14)

冷卻器修正后的有效總阻力系數(shù)如下

(15)

加熱器修正后的有效總阻力系數(shù)如下

(16)

以CFD模型計(jì)算為基準(zhǔn),各種修正方法所得加熱器與冷卻器的波動(dòng)壓降以及時(shí)均聲功的誤差見圖17。

(a)冷卻器波動(dòng)壓降誤差

由圖17可知,經(jīng)修正因子修正后的有效阻力系數(shù)可使Sage一維模型與CFD二維模型計(jì)算的換熱器兩端壓差計(jì)算差異縮小到20%以內(nèi),即Sage一維模型的阻力特性的精度有所提高。換熱器聲功損失修正精度長加熱器比短冷卻器明顯差,這是因?yàn)殚L加熱器的分布特性更難用端部集總壓差簡單修正得到。如表5所示,整體性能對(duì)比表明,Sage一維模型修正后的結(jié)果明顯向CFD二維模型計(jì)算結(jié)果靠近。證明通過有效壓降與有效速度獲得的總阻力系數(shù)可用于修正Sage一維模型的阻力特性,可有效提高工程設(shè)計(jì)精度。

表5 修正前后Sage模型與CFD模型參數(shù)對(duì)比

基于以上分析,可以建立一套提升實(shí)際復(fù)雜結(jié)構(gòu)發(fā)動(dòng)機(jī)的分析效率與精度的方法。首先利用CFD仿真對(duì)真實(shí)換熱器與回?zé)崞鹘M成的核心單元結(jié)構(gòu)進(jìn)行給定聲場(chǎng)下的流動(dòng)阻力特性分析,由于流阻分析的計(jì)算時(shí)間穩(wěn)定較快,可以不需要等待系統(tǒng)熱容與時(shí)均換熱穩(wěn)定即可進(jìn)行分析,因而即使對(duì)于復(fù)雜三維模型的分析計(jì)算同樣具有工程可行性。所得到的換熱器阻力特征可用于修正Sage一維模型中的對(duì)應(yīng)阻力系數(shù)。Sage一維模型中換熱特性的計(jì)算已經(jīng)可以較為準(zhǔn)確地捕捉換熱分布與集中特性,因而無需依賴三維CFD仿真模型的熱穩(wěn)定計(jì)算。時(shí)均換熱特性的CFD計(jì)算由于對(duì)時(shí)間尺度的高要求,以及較大系統(tǒng)熱容穩(wěn)定時(shí)間尺度之間存在的顯著矛盾,使得完全依賴三維CFD仿真進(jìn)行換熱器完整特性分析的計(jì)算時(shí)間過長而難具工程指導(dǎo)意義。

4 結(jié) 論

本文基于10 kW熱管耦合斯特林發(fā)電機(jī)整機(jī)設(shè)計(jì),進(jìn)行發(fā)動(dòng)機(jī)核心部件的拓?fù)涞刃?分別采用Sage一維模型與CFD二維模型兩種計(jì)算方法,對(duì)核心單元的換熱與阻力特性進(jìn)行模擬研究,通過比較兩種模型的計(jì)算結(jié)果,得出以下結(jié)論。

(1)對(duì)換熱器內(nèi)的局部換熱特征以及集總換熱特征計(jì)算結(jié)果對(duì)比,Sage一維模型與CFD二維模型具有相近結(jié)果,Sage一維模型的分析結(jié)果具有工程指導(dǎo)意義。換熱器換熱特性主要由聲場(chǎng)驅(qū)動(dòng),這可能是換熱器換熱特征分析中,對(duì)本構(gòu)方程的準(zhǔn)確性要求相對(duì)較低的原因。

(2)對(duì)換熱器內(nèi)阻力特性的計(jì)算結(jié)果對(duì)比表明,Sage一維模型與CFD二維模型計(jì)算結(jié)果相差較大,換熱器射流導(dǎo)致的回?zé)崞鹘缑骘@著損失是產(chǎn)生這一差異的主要原因。因此,為提高Sage一維計(jì)算結(jié)果的精度,需對(duì)利用CFD對(duì)實(shí)際換熱器與回?zé)崞鹘缑媪髯柽M(jìn)行仿真計(jì)算,提取阻力特性后在Sage一維模型中進(jìn)行修正,即可較好地提高Sage一維模型對(duì)實(shí)際系統(tǒng)的預(yù)測(cè)準(zhǔn)確度。

(3)進(jìn)一步明確了熱聲熱機(jī)換熱器特性的顯著空間分布特征以及工況依賴特征,這些特征在熱聲換熱器中具有相似性和普遍性。但是,對(duì)于特征的量化又因其空間與工況敏感性而難以像穩(wěn)態(tài)流動(dòng)那樣獲得通用的本構(gòu)方程,因而熱聲熱機(jī)換熱器需要針對(duì)具體設(shè)計(jì)進(jìn)行特性分析。