基于體積力的海水自流循環(huán)系統(tǒng)流速確定方法

陳康，姚志崇，周恩東，夏彬

1 深海載人裝備國(guó)家重點(diǎn)實(shí)驗(yàn)室，江蘇無錫214082

2 中國(guó)船舶科學(xué)研究中心，江蘇無錫214082

0 引 言

在現(xiàn)代高速艦船和核動(dòng)力潛艇中，可以采用自流循環(huán)系統(tǒng)為換熱器提供冷卻水［1］。艇體在巡航工況下時(shí)，冷卻系統(tǒng)海水泵不工作而保持隨動(dòng)，利用艇體航行時(shí)與海水相對(duì)運(yùn)動(dòng)產(chǎn)生的動(dòng)水頭來抽吸海水，使海水能克服管路阻力自行流過換熱器，從而有效降低海水泵冷卻系統(tǒng)的負(fù)荷。艇體在高速工況下，需簡(jiǎn)化系統(tǒng)布置，無需另外設(shè)置海水泵；在低速工況下，采用海水泵可輔助實(shí)現(xiàn)海水循環(huán)功能，可降低海水循環(huán)系統(tǒng)對(duì)海水泵功率的需求，從而降低艇體的能耗，提高系統(tǒng)可靠性［2］。

自流循環(huán)系統(tǒng)自流能力的強(qiáng)弱決定了其能否滿足換熱需求，具有良好的自流循環(huán)能力也是設(shè)計(jì)自流循環(huán)系統(tǒng)的根本要求。因此，準(zhǔn)確評(píng)估自流循環(huán)能力（自流循環(huán)系統(tǒng)海水流量或者流速）是設(shè)計(jì)自流循環(huán)系統(tǒng)的關(guān)鍵［3］。由于船身效應(yīng)會(huì)產(chǎn)生壓差，艦船自流循環(huán)系統(tǒng)進(jìn)、出口缺乏有效的工程估算方法，因此無法確定艦船自流循環(huán)系統(tǒng)的自流能力。目前，艦船自流循環(huán)系統(tǒng)的自流能力一般采用CFD 方法進(jìn)行計(jì)算。

裴金亮［4］在運(yùn)用CFD 方法對(duì)自流循環(huán)系統(tǒng)及船體進(jìn)行整體模擬時(shí)，由于自流循環(huán)系統(tǒng)中的換熱器、海水泵、閥門等設(shè)備含有非常細(xì)小的部件，該部件與船體尺度相差較大，從而導(dǎo)致網(wǎng)格劃分存在一定的困難。同時(shí)，模擬自流循環(huán)系統(tǒng)的流場(chǎng)需要布置大量的網(wǎng)格，計(jì)算量很大?？梢?，利用CFD 方法進(jìn)行流量評(píng)估是自流循環(huán)系統(tǒng)設(shè)計(jì)過程中需要考慮的關(guān)鍵問題，由此采用用體積力法等效代替自流循環(huán)系統(tǒng)中部件的方法進(jìn)行模擬具有優(yōu)勢(shì)。

采用力場(chǎng)模擬實(shí)體的體積力法在研究螺旋槳和船體的相互作用時(shí)運(yùn)用較多，其還可以代替真實(shí)的螺旋槳參與仿真計(jì)算。Xu 等［5］和Paterson等［6］運(yùn)用體積力法對(duì)船的伴流場(chǎng)進(jìn)行了計(jì)算。傅慧萍等［7］采用螺旋槳體積力模型，實(shí)現(xiàn)了模型尺度下實(shí)船自航點(diǎn)的全粘帶自由面計(jì)算。吳召華等［8］運(yùn)用一種便捷的描述型體積力模型代替實(shí)體螺旋槳的作用，預(yù)報(bào)了船舶做單平面定常回轉(zhuǎn)運(yùn)動(dòng)時(shí)的受力和力矩。力場(chǎng)模擬方法的計(jì)算量比船槳整體計(jì)算的小，還可以分析槳對(duì)船體的影響，為螺旋槳的設(shè)計(jì)和開發(fā)提供指導(dǎo)。參照體積力法代替螺旋槳的方法，可以將體積力法運(yùn)用到自流循環(huán)系統(tǒng)的數(shù)值模擬中。

海水自流循環(huán)系統(tǒng)的管路阻力、關(guān)鍵部件（換熱器、海水泵、閥門等）阻力都可以通過理論計(jì)算、仿真或?qū)嶒?yàn)得到，由于進(jìn)、出口形狀復(fù)雜，其與艇體的相互作用無法通過理論計(jì)算得到，故預(yù)報(bào)時(shí)只能通過仿真進(jìn)行模擬［9］。為此，本文將仿照體積力等效代替螺旋槳的方法，采用體積力法等效代替換熱器、海水泵和閥門等主要部件，對(duì)自流循環(huán)系統(tǒng)的流場(chǎng)進(jìn)行數(shù)值模擬，建立一種可以確定船舶系統(tǒng)流速和流量的方法，從而為預(yù)報(bào)實(shí)艇海水自流循環(huán)系統(tǒng)的流速和流量提供參考。

1 體積力模擬原理

圖1 所示為某艇的海水自流循環(huán)系統(tǒng)組成示意圖。自流循環(huán)系統(tǒng)一般包括進(jìn)水口部分、換熱器部分、海水泵、中間管路和出水口部分［2］。艇體在水下航行時(shí)，利用航行產(chǎn)生的進(jìn)、出口壓差來驅(qū)動(dòng)海水進(jìn)入自流循環(huán)系統(tǒng)，通過換熱器完成熱量的交換。當(dāng)艇體以低速巡航時(shí)，自流動(dòng)力不足，需要開啟海水泵抽吸海水。自流循環(huán)系統(tǒng)的水頭平衡方程為

式中：vA，vB分別為進(jìn)水口和出水口斷面的平均速度；v0，v1，v2……為自流循環(huán)系統(tǒng)各斷面處的平均速度；PA，PB分別為進(jìn)水口和出水口斷面壓力；ρ為海水密度；g為重力加速度；λ為管路沿程阻力系數(shù)；L為管道長(zhǎng)度；D為管道直徑；ξ進(jìn)為管路進(jìn)口阻力系數(shù)；ξ出為管路出口阻力系數(shù)；ξ換為換熱器阻力系數(shù)；ξ泵為海水泵阻力系數(shù)。

圖1 自流循環(huán)系統(tǒng)組成示意圖Fig.1 Composition of the artesian circulating system

式（1）左邊是艦船以一定速度航行時(shí)，所能提供的總水頭，右邊是管路及各設(shè)備部件的水頭損失，其中換熱器和海水泵的水頭損失為

若已知vA，vB，PA，PB以及各部件的損失系數(shù)，就可以通過式（1）直接求解出自流循環(huán)系統(tǒng)的流速。體積力模擬示意圖如圖2 所示。船體表面為復(fù)雜的三維曲面，自流循環(huán)水系統(tǒng)進(jìn)、出口處的壓力分布沒有準(zhǔn)確的理論求解方法，只能通過仿真計(jì)算方法對(duì)自流循環(huán)系統(tǒng)內(nèi)流場(chǎng)和船體外流場(chǎng)整體進(jìn)行求解來獲取。自流循環(huán)系統(tǒng)內(nèi)各設(shè)備部件可借助體積力來模擬其流動(dòng)損失效應(yīng)，從而避免直接對(duì)各設(shè)備部件進(jìn)行數(shù)值模擬。

圖2 體積力模擬示意圖Fig.2 Body force simulation

為簡(jiǎn)化說明，系統(tǒng)中的部件只選取換熱器。對(duì)計(jì)算域中換熱器部分施加等效換熱器水頭損失的體積力。單位體積的體積力大小為fb，并使

式中：Fb為施加體積力的總值；dV為單位體積。體積力的流動(dòng)控制方程為

在實(shí)際工程中，需要求解的是自流循環(huán)系統(tǒng)的流速，而等效體積力大小是無法用式（2）和式（3）確定的。使用體積力等效代替阻力部件，在對(duì)參考泵選型時(shí)，可以通過2 種特性曲線求交點(diǎn)的方法求得海水自流循環(huán)系統(tǒng)的流速點(diǎn)。水泵的性能曲線如圖3 所示。在水泵的性能曲線上，任意流量點(diǎn)都可以找到與之相對(duì)應(yīng)的揚(yáng)程，該組流量與揚(yáng)程參數(shù)即為工況。管路的特性曲線是管路的水頭損失隨管路流量變化的曲線，因此，管路的特性曲線與水泵性能曲線的交點(diǎn)即為運(yùn)行工況點(diǎn)（即圖3 中的M 點(diǎn)）。

根據(jù)上述方法，通過求取所替代部件的性能曲線與系統(tǒng)特性曲線交點(diǎn)的方法，可以確定實(shí)艇在不同航速下海水自流循環(huán)系統(tǒng)的工況點(diǎn)。系統(tǒng)特性曲線是通過體積力代替主要部件后再與艇體聯(lián)合計(jì)算而求得，部件性能曲線是通過對(duì)所替代部件進(jìn)行單獨(dú)計(jì)算或仿真而得到。

圖3 水泵的性能曲線交點(diǎn)圖Fig.3 Intersection of the pump performance curves

由于用體積力法等效替代阻力部件的方法可以忽略不必要的流動(dòng)細(xì)節(jié)，故在將海水自流循環(huán)系統(tǒng)與大尺寸的艇體聯(lián)合計(jì)算時(shí)，采用體積力法可對(duì)阻力部件進(jìn)行單獨(dú)計(jì)算。同時(shí)，在保證計(jì)算精準(zhǔn)度和有效性的前提下，體積力法還能大幅減少計(jì)算網(wǎng)格，提高計(jì)算效率。

2 自流循環(huán)系統(tǒng)流速確定過程

借助體積力等效代替換熱器等部件來確定自流循環(huán)系統(tǒng)流速的方法，過程如下：

第1 步，在艦船特定航速下，求解自流循環(huán)系統(tǒng)的內(nèi)流場(chǎng)及船體外流場(chǎng)，其中換熱器所在區(qū)域全部設(shè)置為流體域，并分別施加一系列大小不同的體積力fb。計(jì)算求解后，提取換熱器流體域兩端截面的流量和壓力，繪制出流量壓差曲線。

第2 步，由式（2）設(shè)定流速值，求得換熱器的水頭損失，也即海水流經(jīng)換熱器后的壓差。將流速換算為流量，同樣繪制出流量壓差曲線。

第3 步，采用類似于圖3 的方法，求解前2 步繪制的流量壓差曲線的交點(diǎn)，該點(diǎn)即為某航速下自流循環(huán)系統(tǒng)的自流點(diǎn)。第1 步求出的流量壓差曲線相當(dāng)于管路特性曲線，第2 步求出的流量壓差曲線相當(dāng)于水泵性能曲線。讀取自流點(diǎn)對(duì)應(yīng)的流量，就可求出相應(yīng)的流速，也即自流循環(huán)系統(tǒng)的自流流速。

需要進(jìn)一步說明的是，式（2）中換熱器的流動(dòng)損失系數(shù)ξ換需根據(jù)經(jīng)驗(yàn)或?qū)嶒?yàn)確定。而根據(jù)艦船自流循環(huán)系統(tǒng)設(shè)計(jì)的需求，換熱器一般無法使用標(biāo)準(zhǔn)換熱器。在設(shè)計(jì)非標(biāo)準(zhǔn)換熱器時(shí)，其流動(dòng)損失系數(shù)未知，這時(shí)，可以采用CFD 方法確定。使用該方法確定流動(dòng)損失系數(shù)，只需要單獨(dú)模擬換熱器，最重要的是，還能與自流循環(huán)系統(tǒng)及船體設(shè)計(jì)同步進(jìn)行，且不影響系統(tǒng)多方案設(shè)計(jì)的對(duì)比，如自流循環(huán)系統(tǒng)進(jìn)、出口流道的設(shè)計(jì)以及位置選擇等。在海水自流循環(huán)系統(tǒng)上采用體積力法有2個(gè)優(yōu)勢(shì)：一是換熱器采用CFD 模擬法可用于自流循環(huán)系統(tǒng)多方案設(shè)計(jì)中，而不是每設(shè)計(jì)一種換熱器方案都要對(duì)整艇進(jìn)行數(shù)值模擬；二是換熱器和自流循環(huán)系統(tǒng)的設(shè)計(jì)可同步進(jìn)行，從而縮短設(shè)計(jì)周期。

其他設(shè)備部件也可以采用體積力法進(jìn)行等效代替。為簡(jiǎn)化說明，后續(xù)實(shí)例中的部件只考慮換熱器，需要加入其他部件時(shí)只需要將部件的特性曲線并入到總的部件性能曲線中即可。

3 CFD 數(shù)值模擬方法

3.1 網(wǎng)格劃分

整個(gè)自流循環(huán)系統(tǒng)的研究對(duì)象不僅包括系統(tǒng)內(nèi)部的內(nèi)流場(chǎng)，還包括艇體外部的外流場(chǎng)。為了模擬整個(gè)自流循環(huán)系統(tǒng)的內(nèi)部流動(dòng)特性，參加計(jì)算的部分包括外流場(chǎng)、進(jìn)出水管道以及換熱器。海水自流循環(huán)系統(tǒng)模擬的內(nèi)、外流場(chǎng)的流動(dòng)速度與幾何尺寸相差較大，導(dǎo)致其流動(dòng)特性也存在很大的不同。使用GAMBIT 軟件分塊畫網(wǎng)格，然后再組合到一塊進(jìn)行內(nèi)、外流場(chǎng)聯(lián)合計(jì)算。實(shí)際選用圓柱形流體計(jì)算域，流場(chǎng)半徑為艇體型深的12倍；在長(zhǎng)度方向上，計(jì)算域向艇體艏部延伸艇長(zhǎng)的2 倍，向艇體艉部延伸艇長(zhǎng)的3 倍。自流循環(huán)系統(tǒng)流體域的網(wǎng)格劃分如圖4 所示。

圖4 模型網(wǎng)格劃分Fig.4 Model mesh division

在使用體積力法等效代替換熱器進(jìn)行數(shù)值模擬時(shí)，代入的是體積力的總值而不是系統(tǒng)的流場(chǎng)分布。因此，體積力的分布方式不影響自流循環(huán)系統(tǒng)流速點(diǎn)的求取，實(shí)際將體積力均勻分布到圓柱形流體域上。體積力法計(jì)算示意圖如圖5所示。圖中：A，B 分別為海水進(jìn)、出口截面；C，D 為施加體積力的圓柱形區(qū)域邊界的截面。體積力代替部分網(wǎng)格數(shù)為1 萬，整個(gè)流場(chǎng)的網(wǎng)格總數(shù)為800 萬。

為驗(yàn)證體積力法等效代替部件方法的準(zhǔn)確性，將體積力法與傳統(tǒng)計(jì)算方法的結(jié)果進(jìn)行對(duì)比。將換熱器實(shí)體劃分網(wǎng)格后加到自流循環(huán)系統(tǒng)流道中，然后再進(jìn)行整艇的聯(lián)合模擬計(jì)算。實(shí)體聯(lián)合計(jì)算網(wǎng)格如圖6 所示，最后采用的換熱器網(wǎng)格數(shù)為1 036 萬，整個(gè)流場(chǎng)網(wǎng)格總數(shù)為1 832 萬。為了在模擬艇體外流場(chǎng)的同時(shí)兼顧自流循環(huán)系統(tǒng)的內(nèi)流場(chǎng)，于計(jì)算后再劃分網(wǎng)格，使換熱器的壁面y+≈4，艇體艉部和自流循環(huán)系統(tǒng)進(jìn)、出口的壁面y+≈102，艇體和外流場(chǎng)其他區(qū)域的壁面y+≈275。這里，用y+表征網(wǎng)格密度。

圖5 體積力法計(jì)算示意圖Fig.5 Calculation diagram of body force method

圖6 實(shí)體聯(lián)合計(jì)算示意圖Fig.6 Combined calculation diagram of component and hull

自流循環(huán)系統(tǒng)使用的換熱器以管殼式換熱器為例。劃分管殼式換熱器網(wǎng)格的難點(diǎn)在于，換熱器的換熱管太過細(xì)長(zhǎng)，與換熱器整體尺度相比較小。由于系統(tǒng)部件的尺寸變化大，為保證計(jì)算的精準(zhǔn)度，需要在換熱管的內(nèi)壁面劃分大量網(wǎng)格。在劃分換熱器網(wǎng)格時(shí)，為了能更好地對(duì)管壁處的流場(chǎng)進(jìn)行模擬，需要對(duì)壁面附近的網(wǎng)格進(jìn)行加密處理。換熱器網(wǎng)格如圖7 所示。

圖7 換熱器網(wǎng)格圖Fig.7 Heat exchanger grid

3.2 湍流模型

仿真計(jì)算涉及艇體的外部流場(chǎng)和內(nèi)部海水管路流場(chǎng)。在外流場(chǎng)和內(nèi)流場(chǎng)的混合計(jì)算過程中，需采用合適的湍流模型。FLUENT 軟件提供的兩方程湍流模型主要包括標(biāo)準(zhǔn)k-ε模型、RNG k-ε模型、Realizable k-ε模型、標(biāo)準(zhǔn)k-ω模型和SST k-ω模型［10］。張楠等［11］結(jié)合這幾種湍流模型，將計(jì)算得到的潛艇阻力和尾流場(chǎng)與模型試驗(yàn)結(jié)果進(jìn)行了對(duì)比分析，驗(yàn)證了數(shù)值方法的可靠性。張志榮等［12］針對(duì)這幾種不同湍流模型的粘性繞流場(chǎng)進(jìn)行了數(shù)值模擬，將計(jì)算得出的阻力系數(shù)、槳盤面處流場(chǎng)與試驗(yàn)值進(jìn)行了比較。結(jié)果表明，在復(fù)雜船型粘性流場(chǎng)計(jì)算的湍流模式中，和標(biāo)準(zhǔn)k-ω模型相比，使用SST k-ω模型具有更高的精準(zhǔn)度和可靠性，其綜合了k-ω模型在近壁面區(qū)域和遠(yuǎn)場(chǎng)區(qū)域計(jì)算的優(yōu)點(diǎn)，同時(shí)增加了橫向耗散導(dǎo)數(shù)項(xiàng)，在湍流粘度定義中考慮了湍流剪切應(yīng)力的運(yùn)輸過程，更加適用于廣域流場(chǎng)。因此，本文實(shí)際選用的是SST k-ω模型進(jìn)行計(jì)算。

3.3 數(shù)值計(jì)算方法

數(shù)值模擬求解三維定常不可壓縮RANS 方程，采用有限體積法離散動(dòng)量方程，擴(kuò)散項(xiàng)使用中心差分格式，對(duì)流項(xiàng)采用二階迎風(fēng)格式，壓力速度耦合采用SIMPLE 算法，離散得到的代數(shù)方程組用Gauss-Seidel 法迭代求解，并使用代數(shù)多重網(wǎng)格法加速收斂。

3.4 體積力加載

將體積力均勻施加到以換熱器直徑、長(zhǎng)度所定義的圓柱形流體域內(nèi)，來實(shí)現(xiàn)實(shí)際部件對(duì)海水的阻力影響。按照表1 所示數(shù)值施加體積力，選定某特定航速后，通過施加對(duì)應(yīng)的體積力，可以模擬得到該工況下系統(tǒng)的流速和流量。體積力施加方式如圖2 所示，體積力施加方向與海水流速方向相反。

表1 體積力施加值Table 1 Values of the applied body force at 2，4，6 knots

4 自流循環(huán)流速點(diǎn)的確定

4.1 系統(tǒng)特性曲線

圖5 中，C，D 截面的壓差隨流速變化的曲線即為體積力代替部件后的系統(tǒng)特性曲線。選定航速后，施加表1 中航速對(duì)應(yīng)的體積力進(jìn)行模擬計(jì)算，3 種航速下不同體積力對(duì)應(yīng)的海水管道流速如圖8 所示。由圖8 可以看出，在一定航速下，隨著體積力的增加，自流循環(huán)的流速是下降的；在相同體積力下，隨著航速的增加，自流循環(huán)的流速和流量是不斷增加的。系統(tǒng)特性曲線如圖9 所示，從圖9 中可以看出，在相同海水管道流速和流量下，隨著航速的增加，C，D 截面的壓差逐漸增大，符合實(shí)際趨勢(shì)。

圖8 在2，4，6 kn 航速下不同體積力對(duì)應(yīng)海水管道流速Fig.8 Velocities of sea water pipe at 2，4，6 knots in different body forces

圖9 在2，4，6 kn 航速下不同流速對(duì)應(yīng)C，D 截面的壓差Fig.9 Differential pressure between section C and D at 2，4，6 knots

4.2 部件性能曲線

自流循環(huán)系統(tǒng)的部件性能曲線為不同管道流速下流體域等效代替部件的壓差。此處流體域等效代替的部件為管殼式換熱器，換熱器的管程壓降可以通過工程理論［13］或單獨(dú)仿真計(jì)算得到。在需要添加其他阻力部件時(shí)，可通過相同的方法計(jì)算或仿真以后，再并入到總的部件性能曲線當(dāng)中。自流循環(huán)系統(tǒng)的阻力計(jì)算如表2 所示。

表2 部件壓差計(jì)算Table 2 Differential pressure of the component

4.3 求取曲線交點(diǎn)

將系統(tǒng)特性曲線和部件性能曲線繪制到同一張圖中，曲線的交點(diǎn)即為海水自流循環(huán)系統(tǒng)在特定航速下的流速點(diǎn)。圖中兩條曲線存在交點(diǎn)說明在該航速下可以實(shí)現(xiàn)海水的自流循環(huán)，若沒有交點(diǎn)，則說明此種工況下的自流循環(huán)系統(tǒng)流速過低，不能滿足冷卻水量的需求，需要開啟海水泵主動(dòng)循環(huán)系統(tǒng)，以增加吸水揚(yáng)程。部件性能曲線與系統(tǒng)特性曲線的交點(diǎn)圖如圖10 所示。

圖10 系統(tǒng)特性曲線與部件性能曲線交點(diǎn)圖Fig.10 Intersection of system characteristic curves and component performance curve

5 聯(lián)合計(jì)算驗(yàn)證

為了驗(yàn)證體積力法的準(zhǔn)確性，將聯(lián)合計(jì)算法與體積力法計(jì)算得到的自流循環(huán)工況點(diǎn)進(jìn)行了比較。單獨(dú)劃分換熱器的網(wǎng)格后，將換熱器實(shí)體加入到自流循環(huán)管路中與艇體一起進(jìn)行數(shù)值模擬，仿真計(jì)算的參數(shù)設(shè)置與體積力法的參數(shù)設(shè)置相同。選取2，4，6 kn 航速下的工況，計(jì)算得到自流循環(huán)系統(tǒng)的流速和部件兩端的壓差，體積力法與聯(lián)合計(jì)算法的計(jì)算結(jié)果如表3 所示。

將體積力法與實(shí)體聯(lián)合計(jì)算法得到的自流循環(huán)工況點(diǎn)進(jìn)行比較，結(jié)果如圖11 所示。由圖可看出，2 種方法下的流速值與壓差值相差均在5%以內(nèi)。對(duì)比2 種方法的網(wǎng)格總數(shù)，發(fā)現(xiàn)體積力法具有顯著優(yōu)勢(shì)。

表3 自流循環(huán)工況點(diǎn)比較Table 3 Comparison of the operating point

圖11 體積力法與聯(lián)合計(jì)算法工況點(diǎn)比較Fig.11 Comparison of operating point between body force method and combined calculation

6 結(jié) 語

本文基于體積力法建立了確定自流循環(huán)系統(tǒng)流速的分析方法，運(yùn)用體積力法對(duì)自流循環(huán)系統(tǒng)進(jìn)行數(shù)值模擬，類似于確定水泵工況點(diǎn)的方法，通過求取系統(tǒng)特性曲線與部件性能曲線交點(diǎn)的方法得到海水自循環(huán)系統(tǒng)的流速點(diǎn)。通過比較實(shí)體聯(lián)合計(jì)算法和體積力法這2 種方法下的自流循環(huán)工況點(diǎn)，表明體積力模擬方法可以確定海水自循環(huán)系統(tǒng)的流速，而且計(jì)算量小。后續(xù)若需添加其他阻力部件，可以直接將對(duì)應(yīng)的部件性能曲線并入到總的部件性能曲線當(dāng)中，通過求交點(diǎn)的方法來得出實(shí)際海水流量。有關(guān)仿真結(jié)果的準(zhǔn)確性還有待后續(xù)的試驗(yàn)驗(yàn)證。