基于CFD的雙向流道泵站壓力脈動混沌特性研究

肖忠明，顏紅勤，蔣紅櫻，成 立，劉志泉

（1.丹陽市九曲河樞紐管理處，江蘇鎮(zhèn)江 212000；2.江蘇省水利工程科技咨詢股份有限公司，南京 210029；3.揚(yáng)州大學(xué)水利科學(xué)與工程學(xué)院，江蘇揚(yáng)州 225009）

0 引 言

在我國長江沿岸，為滿足灌溉、排澇等要求，建有許多雙向流道泵站。該類泵站的特點(diǎn)是流量大、揚(yáng)程低［1］，同時受長江水位變化影響，易偏離設(shè)計(jì)工況運(yùn)行［2，3］。水泵長期偏離設(shè)計(jì)工況運(yùn)行，不僅會誘發(fā)不良壓力脈動［3，4］，導(dǎo)致機(jī)組異常振動，而且降低機(jī)組運(yùn)行效率，影響泵站安全穩(wěn)定性，引發(fā)泵站故障［5-7］。

針對軸流泵機(jī)組壓力脈動特性，國內(nèi)外許多學(xué)者進(jìn)行了研究［8-12］。王福軍、施衛(wèi)東等［8，9］采用數(shù)值模擬的方法研究了軸流泵在不同工況下的壓力脈動特性，得出壓力脈動幅值在葉輪進(jìn)、出口隨半徑增大而增大，導(dǎo)葉中間和出口變化規(guī)律則相反，且偏離最優(yōu)工況越遠(yuǎn)，脈動的相對振幅越大。張德勝等［10］采用試驗(yàn)的方法研究了軸流泵在不同流量、轉(zhuǎn)速下的壓力脈動特性，揭示了軸流泵內(nèi)部不同位置處壓力脈動規(guī)律。Al-Obaidi A R［11］研究了軸流泵在不同葉片角度下壓力脈動的時域和頻域特性，結(jié)果表明葉輪葉片角度和非定常流動對壓力脈動幅值影響較大。目前對軸流泵壓力脈動特性的研究分析方法主要是FFT、小波分析等，此類方法對脈動信號周期性明顯的葉輪區(qū)域有較好的適用性。但在遠(yuǎn)離振動源葉輪的區(qū)域，如導(dǎo)葉出口、進(jìn)出水喇叭口、流道等，壓力脈動分布規(guī)律逐漸變得雜亂無序［13，14］，此類方法適用性較差?；煦缋碚撜J(rèn)為在非線性系統(tǒng)中看似隨機(jī)不確定性的行為，不是來源于外部環(huán)境的隨機(jī)因素對系統(tǒng)運(yùn)動的影響，而是系統(tǒng)自發(fā)產(chǎn)生的，即具有隨機(jī)性的運(yùn)動狀態(tài)是由確定性方程（內(nèi)因）直接得到的［15］。許多學(xué)者利用混沌理論對泵和水輪機(jī)等流體機(jī)械進(jìn)行壓力脈動分析和故障預(yù)測［16-22］，梁超等［16］對采集的壓力脈動信號進(jìn)行小波包分解，根據(jù)混沌特征參數(shù)把離心泵汽蝕過程劃分為4 個階段，從不同的角度定量描述動力系統(tǒng)的混沌特性。蘇文濤等［17］研究了水輪機(jī)在不同空化程度下的壓力脈動信號的混沌動力學(xué)特征，認(rèn)為水輪機(jī)壓力脈動信號中存在混沌吸引子。

以某立式軸流泵站為模型，通過CFD 計(jì)算得到設(shè)計(jì)工況下各監(jiān)測點(diǎn)的壓力脈動時域信號，對其進(jìn)行相空間重構(gòu)，利用小數(shù)據(jù)量法求得最大李雅普諾夫指數(shù)，判斷其混沌特性，通過求取其分形標(biāo)度和關(guān)聯(lián)維數(shù)，對壓力脈動信號的混沌特性進(jìn)行分析，并對其他兩種常用工況進(jìn)行驗(yàn)證。

1 數(shù)值計(jì)算模型及計(jì)算方法

1.1 計(jì)算模型

以某沿江雙向流道泵站為計(jì)算模型，如圖1所示，該泵站設(shè)計(jì)流量Qopt=20 m3/s，設(shè)計(jì)揚(yáng)程H=2.61 m，葉輪直徑D=2 550 mm，轉(zhuǎn)速n=150 r/min，葉片數(shù)為3片，導(dǎo)葉片數(shù)為7片，采用“X”型雙層箱涵式進(jìn)出水流道。

圖1 計(jì)算模型示意圖Fig.1 Diagram of calculation model

1.2 網(wǎng)格劃分

對計(jì)算模型進(jìn)行六面體結(jié)構(gòu)化網(wǎng)格劃分，對葉輪、導(dǎo)葉及邊壁進(jìn)行加密處理，劃分網(wǎng)格數(shù)量962.6 萬個，整體網(wǎng)格質(zhì)量控制在0.5以上。葉輪和導(dǎo)葉網(wǎng)格劃分示意圖如圖2所示。

1.3 邊界條件設(shè)置

CFD 計(jì)算采用非定常計(jì)算，計(jì)算模型進(jìn)口設(shè)為質(zhì)量流量進(jìn)口，出口設(shè)為自由出流條件；葉輪與導(dǎo)葉、進(jìn)水流道的交界面采用瞬態(tài)動靜轉(zhuǎn)子模型。文獻(xiàn)［23］的研究結(jié)果表明，SSTk-ω模型對流體機(jī)械的計(jì)算有較好的適用性，因此湍流模型選擇SSTk-ω模型，計(jì)算收斂精度為10-5。葉輪旋轉(zhuǎn)周期T=0.4 s，取葉輪旋轉(zhuǎn)周期的1/90 為1 個時間步長，即時間步長Δt=3.333 3×10-3s，共計(jì)算12個周期。

1.4 測點(diǎn)布置

為監(jiān)測軸流泵內(nèi)部壓力脈動，在進(jìn)水喇叭口、葉輪進(jìn)口、葉輪出口、導(dǎo)葉中部、導(dǎo)葉出口、出水喇叭口截面各布置3個測點(diǎn)，從輪緣至輪轂沿徑向分布，共布置18 個測點(diǎn)，測點(diǎn)布置示意圖如圖2所示。

圖2 網(wǎng)格劃分及測點(diǎn)布置示意圖Fig.2 Diagram of grid division and measuring point arrangement

2 外特性驗(yàn)證

在CFX 中對該泵站在0°葉片安放角下的7 個工況進(jìn)行計(jì)算，流 量 分 別 為0.7Qopt、0.8Qopt、0.9Qopt、1.0Qopt、1.1Qopt、1.2Qopt、1.3Qopt，得到其關(guān)于流量Q和揚(yáng)程H、效率η的外特性曲線，并與試驗(yàn)換算得到的原型泵裝置性能數(shù)據(jù)對比，如圖3所示。

圖3 數(shù)模與試驗(yàn)外特性曲線圖Fig.3 Curves of numerical simulation and test external characteristics

由圖3 可知，數(shù)值模擬得到的外特性曲線比試驗(yàn)外特性曲線整體略低，兩者變化趨勢基本一致。因此，SSTk-ω模型可以較好地預(yù)測軸流泵的內(nèi)部流動，數(shù)值模擬計(jì)算結(jié)果可靠。

3 壓力脈動時域信號分析

通過非定常計(jì)算得到泵站運(yùn)行時的壓力脈動時域信號，選取徑向外緣有代表性的測點(diǎn)P3、P6、P9、P12、P15、P18進(jìn)行分析，測點(diǎn)的壓力脈動時域圖如圖4 所示，為方便觀察，圖中均顯示后4個周期的時域圖。

由圖4（a）可知，進(jìn)水喇叭口處P3點(diǎn)壓力脈動波形較亂，無明顯的周期性，因?yàn)樗鲝牧鞯肋M(jìn)入進(jìn)水喇叭口時，水流速度、方向改變，使幅值呈現(xiàn)波動變化。由圖4（b）、圖4（c）可知，葉輪進(jìn)口P6點(diǎn)周期性規(guī)律明顯，1 個周期內(nèi)存在1 個波峰波谷；葉輪出口P9點(diǎn)壓力脈動波形復(fù)雜，1個周期內(nèi)有多個波峰波谷，幅值變化較大，此處水流受葉輪和導(dǎo)葉動靜干涉作用的影響，流場變化復(fù)雜。由圖4（d）可知，導(dǎo)葉中部P12點(diǎn)壓力脈動變化受葉輪的影響，1 個周期內(nèi)有3 個波形波谷，與葉片數(shù)相對應(yīng)，幅值發(fā)生變化可能與導(dǎo)葉內(nèi)存在回流等不良流動有關(guān)。由圖4（e）、圖4（f）可知，導(dǎo)葉出口P15點(diǎn)、出水喇叭口P18點(diǎn)距葉輪較遠(yuǎn)，此處壓力脈動規(guī)律不明顯，但幅值較高，雙向出水流道內(nèi)流態(tài)較差，對該處壓力脈動的變化產(chǎn)生影響。

圖4 測點(diǎn)壓力脈動時域圖Fig.4 Time-domain diagram of pressure at each measuring point

4 混沌特性分析

4.1 相空間重構(gòu)

混沌時間序列的判定、分析與預(yù)測都是在重構(gòu)的相空間中進(jìn)行的，因此相空間重構(gòu)是混沌特性分析重要的步驟［24］。采用坐標(biāo)延遲法［24，25］將計(jì)算得到的壓力脈動時間序列｛x（t）｝，t=1，2，…，N，根據(jù)時間延遲τ與嵌入維數(shù)m進(jìn)行相空間重構(gòu)，得到一組新序列：

式中：t=1，2，…，M；M=N-（m-1）τ。

時間延遲τ采用互信息量法求得，嵌入維數(shù)m采用G-P 算法求得，詳細(xì)計(jì)算方法參見文獻(xiàn)［24，26，27］。選取有代表性的測點(diǎn)P3，P6，P9，P12，P15，P18點(diǎn)的壓力脈動進(jìn)行相空間重構(gòu)，得到三維相軌跡圖如圖5所示。

由圖5可知，進(jìn)水喇叭口P3點(diǎn)相軌跡圖無明顯規(guī)律，從底部至上部呈現(xiàn)擴(kuò)張規(guī)律。根據(jù)葉輪進(jìn)口P6點(diǎn)壓力脈動時域圖4（b），該點(diǎn)壓力脈動周期性最明顯，相軌跡圖呈收縮狀，表現(xiàn)出混沌特性。葉輪出口P9點(diǎn)相軌跡圖呈扭轉(zhuǎn)的圓環(huán)狀，混沌特性明顯，對應(yīng)該點(diǎn)復(fù)雜的壓力脈動變化。導(dǎo)葉中部P12點(diǎn)相軌跡圖呈現(xiàn)螺旋上升的曲線狀，該點(diǎn)壓力脈動仍具有周期性，但幅值有所變化，使相軌跡圖收縮為長條狀，序列點(diǎn)自底部至上部離散程度增加。導(dǎo)葉出口P15點(diǎn)、出水喇叭口P18點(diǎn)遠(yuǎn)離振動源，加之受流道水流變化的影響，壓力脈動無周期性，相軌跡圖為不規(guī)則的長條狀。

圖5 測點(diǎn)壓力脈動三維相軌跡圖（無量綱）Fig.5 Three-dimensional phase trajectory of pressure pulsation at each measuring point

4.2 最大李雅普諾夫指數(shù)分析

李雅普諾夫指數(shù)（LE）是混沌的表征量之一，可用李雅普諾夫指數(shù)定量分析偏流量工況下的壓力脈動信號，判別和驗(yàn)證系統(tǒng)的混沌程度。采用小數(shù)據(jù)量法［26，27］求解設(shè)計(jì)流量下的P3、P6、P9、P12、P15、P18號點(diǎn)的最大李雅普諾夫指數(shù)。在混沌系統(tǒng)中最大李雅普諾夫指數(shù)λ＞0，表示系統(tǒng)既不會穩(wěn)定在不動點(diǎn)，也不存在穩(wěn)定的周期解，同時也不會發(fā)散，系統(tǒng)進(jìn)入混沌。求得各監(jiān)測點(diǎn)脈動信號的最大李雅普諾夫指數(shù)均大于0，分別為0.000 23、0.000 39、0.000 24、0.000 086、0.001 506、0.002 235，因此該雙向流道泵機(jī)組壓力脈動信號具有混沌特性。

4.3 關(guān)聯(lián)維分析

關(guān)聯(lián)維數(shù)作為分形維數(shù)的一種，能描述時間序列相空間中奇異吸引子的維數(shù)，可由時間序列直接得到［28］。關(guān)聯(lián)維數(shù)度量了重構(gòu)相空間的向量距離隨超球半徑變化的速率，也度量了系統(tǒng)內(nèi)部結(jié)構(gòu)復(fù)雜程度。采用G-P 算法計(jì)算關(guān)聯(lián)維數(shù)，關(guān)聯(lián)積分C（r）定義如下：

式中：r為相空間中超球體半徑；Yi為重構(gòu)相空間中的向量。

作出關(guān)聯(lián)積分曲線lnr～lnC（r），對曲線上的直線區(qū)域擬合得到關(guān)聯(lián)曲線。逐漸增加嵌入維數(shù)m直至其斜率不發(fā)生明顯變化，即呈現(xiàn)飽和趨勢，對此時的斜率進(jìn)行擬合即得到關(guān)聯(lián)維數(shù)D，此時的m為最優(yōu)嵌入維數(shù)。

圖6 為設(shè)計(jì)工況下P3、P6、P9、P12、P15、P18號監(jiān)測點(diǎn)嵌入維數(shù)從1遞增到30的關(guān)聯(lián)積分曲線圖，擬合即得關(guān)聯(lián)維數(shù)D，分別為2.7、2.08、3.3、1.96、2.02、1.96。由圖6 可知，隨著嵌入維數(shù)m的增大，關(guān)聯(lián)維數(shù)D也相應(yīng)增大，直至趨于飽和，近乎不變，這是混沌時間序列的特征，設(shè)計(jì)工況下各監(jiān)測點(diǎn)壓力脈動信號的混沌特性得以驗(yàn)證。葉輪出口P9點(diǎn)處監(jiān)測點(diǎn)的關(guān)聯(lián)維數(shù)最大，出水喇叭口P18點(diǎn)處關(guān)聯(lián)維數(shù)最小。結(jié)合圖4 來看，關(guān)聯(lián)維數(shù)大小與壓力脈動振幅變化大小相對應(yīng)，振幅變化越大，關(guān)聯(lián)維數(shù)越大；振幅變化越小，關(guān)聯(lián)維數(shù)越小。

圖6 各測點(diǎn)關(guān)聯(lián)積分曲線圖Fig.6 Correlation integral curve of each measuring point

4.4 分形標(biāo)度分析

利用去傾擾動分析法（DFA）［26］對分形標(biāo)度進(jìn)行求值。對壓力脈動時間序列｛x（t）｝求平均值xˉ，分別與均值作差，累加得序列y（k）：

將y（k）平均分為l份的新序列，每一份的長度即為n=N/l。通過最小二乘法對新序列進(jìn)行擬合，求得每個子序列的擬合函數(shù)yn（k）后，求取y（k）與yn（k）差值的均方根F（n）：

選取每個監(jiān)測點(diǎn)的后800個數(shù)據(jù)，取100以內(nèi)800的因數(shù)得n=5，8，10，16，20，25，32，40，50，80，100，進(jìn)行l(wèi)gn-lgF（n）曲線圖的繪制如圖7 所示，擬合的對數(shù)曲線的斜率即為相應(yīng)的時間序列對應(yīng)的分形標(biāo)度值。

由圖7 可知，隨著lgn的增大，lgF（n）也隨之增大，擬合曲線的斜率逐漸減小，求得的相應(yīng)的時間序列對應(yīng)的分形標(biāo)度值分別為1.683 8、0.722 1、2.195 1、1.957 2、2.546 4、1.842 0，從數(shù)據(jù)上可以看出分形標(biāo)度值均大于0，可知每個時間序列均為混沌時間序列。監(jiān)測點(diǎn)P6分形標(biāo)度值最接近0，最接近臨界混沌狀態(tài)；監(jiān)測點(diǎn)P15分形標(biāo)度值最大，代表在此情況下混沌特性最強(qiáng)。

圖7 設(shè)計(jì)工況F（n）-lg n曲線圖Fig 7 Diagramlog F（n）-lg n of design condition

5 結(jié) 論

（1）數(shù)模得到的外特性曲線與試驗(yàn)外特性曲線變化趨勢一致，表明計(jì)算結(jié)果可靠。

（2）進(jìn)水喇叭口處水流、速度方向發(fā)生變化，壓力脈動波形較亂，無明顯的周期性；葉輪進(jìn)口受葉輪轉(zhuǎn)動影響，壓力脈動周期性規(guī)律明顯；葉輪出口受動靜干涉作用，壓力脈動波形復(fù)雜；導(dǎo)葉出口、出水喇叭口距葉輪較遠(yuǎn)，壓力脈動規(guī)律不明顯，但幅值較高。

（3）各測點(diǎn)最大李雅普諾夫指數(shù)均大于0，表明雙向流道泵站壓力脈動時域信號具有混沌特性，可用混沌理論進(jìn)行分析研究。

（4）隨著嵌入維數(shù)的增加，關(guān)聯(lián)維數(shù)呈現(xiàn)飽和的趨勢，進(jìn)一步證明壓力脈動時域信號具有混沌特性，計(jì)算得到的時間序列均為混沌時間序列。P9測點(diǎn)振幅變化最大，關(guān)聯(lián)維數(shù)最大；P18測點(diǎn)振幅變化最小，關(guān)聯(lián)維數(shù)最??；各監(jiān)測點(diǎn)關(guān)聯(lián)維數(shù)大小與壓力脈動振幅變化大小相對應(yīng)。

（5）各時間序列對應(yīng)的分形標(biāo)度值均大于0，監(jiān)測點(diǎn)P6分形標(biāo)度值最接近0，表明該點(diǎn)最接近臨界混沌狀態(tài)；監(jiān)測點(diǎn)P15分形標(biāo)度值最大，表明該點(diǎn)下混沌特性最強(qiáng)?；煦缋碚摓榉治霰谜镜膲毫γ}動提供了一種新方法。