離心泵葉輪內非定常流動的動態模態分解分析

陳學炳， 張人會,2， 蔣利杰， 郭廣強

(1. 蘭州理工大學 能源與動力工程學院，蘭州 730050； 2. 蘭州理工大學 甘肅省流體機械及系統重點實驗室，蘭州 730050)

離心泵作為一種通用機械，被廣泛應用于各個領域，對國民經濟具有重要作用。復雜的非定常流動是造成離心泵內壓力脈動、振動、噪聲等不穩定流動現象的主要原因，這些不穩定流動往往會造成泵內較大的水力損失，甚至危及泵的安全穩定運行[1-4]。因此，對離心泵內非定常流場特性進行研究具有重要意義。目前，對于離心泵內復雜的非穩態流動，國內外研究人員已開展了大量的工作[5-10]。高波等[11]基于數值模擬方法對離心泵葉輪出口沿軸向變化的圓周面流域內各節點的壓力分布特性進行分析，發現在不同工況下葉輪出口圓周方向的壓力系數均在隔舌附近產生極值。Pedersen等[12]和Byskov等[13]分別采用試驗與數值模擬的方法對離心泵葉輪內流動進行研究，并將數值模擬結果與PIV試驗測試結果進行對比分析。Zhang等[14]采用延遲分離渦模型(delay detached eddy simulation, DDES)對離心泵內流場進行模擬，詳細分析了葉片掃掠蝸殼隔舌時的非定常尾流及其在不同時刻的演變規律。郭榮等[15]對射流式離心泵過流部件誘發的流動噪聲和流激噪聲特性進行分析，發現葉輪和導葉的動靜干涉是誘發射流式離心泵內場噪聲的重要因素之一。離心泵內部復雜流場具有多尺度的流動特性，目前對其復雜瞬態流動的研究大多局限于對其壓力脈動和泵內部整體流場的測量及分析，但對于離心泵葉輪內非定常流動尤其是相對流動中不同頻率所對應的流場結構特征及其相應的物理機理方面的研究較少，若可以將其復雜流場分解為若干不同尺度的流動結構，將對葉輪內流動機理研究提供一定的參考。

動態模態分解(dynamic modal decomposition, DMD)是一種從非定常試驗測量或數值模擬流動中提取動力學信息的數據驅動算法[16-20]。Schmid[21-22]首次提出DMD方法，并將其用于對高維復雜流場進行低維近似，同時識別出脈動主頻及其對應的空間結構。Rowley等[23]采用DMD方法對橫向射流中的非線性流動進行分析，進而提取出不同頻率的流場模態特征。洪樹立等[24]采用DMD方法將包含復雜時空信息的葉柵分離流流場進行解耦，獲得了流場的主要頻率及與之對應的空間結構。張人會等[25]對螺旋軸流式油氣混輸泵導葉內流場進行模態分析，系統研究了復雜流動動態模態的流場結構及其對應的物理機理。吳亞東等[26]采用DMD方法對軸流式壓氣機的旋轉不穩定性進行分析，將其流場按照頻率分解為成對出現的不同模態，且各階模態從頻率和空間周期上與不同的空間周向模態相對應。

目前的研究表明，DMD方法在流場分解方面取得了一定的進展，但基于DMD方法分析葉輪內相對流動的非定常特性仍需進一步研究。本文基于離心泵內流場的大渦模擬(large eddy simulation, LES)數值計算結果，采用DMD方法對其葉輪內部復雜非定常流場進行動態模態分解，分析其瞬態相對速度場在不同頻率下的流場結構及其隨時間演化的動態特性，以研究葉輪內復雜流場中包含的流動特征及其對應的物理機制。

1 數值計算

1.1 數值計算模型

本研究以一臺中比轉速單級單吸離心泵為研究對象，其基本設計參數為：流量Q=200 m3/h，揚程H=20 m，轉速n=1 450 r/min，效率η=81%，比轉速ns=132，葉輪進口直徑D1=150 mm，葉輪外徑D2=260 mm，葉輪出口寬度b=33 mm，葉片數Z=6，蝸殼基圓直徑D3=262 mm，蝸殼出口直徑D4=125 mm。數值計算域包括進口段、葉輪、蝸殼及出口段。采用六面體結構化網格對整個數值計算域進行離散，對近壁面等流動參數變化較大的區域進行局部加密，第一層近壁面網格距離壁面0.012 mm，葉片表面平均y+值約為1.5，其他區域平均y+值約為5。采用定常模擬方法對設計工況下不同網格數量的4種方案進行數值計算，如表1所示，以驗證網格無關性。隨著網格數量的增加泵揚程逐漸趨于穩定，綜合考慮計算資源及數值計算精度，最后確定網格總數為650萬，數值計算網格如圖1所示。

表1 網格無關性測試Tab.1 Grid independence test

圖1 數值計算網格Fig.1 Grid of numerical calculation

1.2 模型邊界條件設置

采用ANSYS Fluent 16.0對離心泵內部復雜流動進行數值模擬。定常數值計算采用SSTk-ω湍流模型，速度進口及自由出流邊界條件，壁面滿足無滑移固壁條件。為更加精確地描述離心泵內部復雜流動，非定常計算采用大渦模擬進行數值求解，亞格子模型選擇Smagorinsky-Lilly模型，動靜交界面采用滑移網格處理。采用SIMPLEC算法耦合壓力及速度場，收斂精度設置為10-5，動量方程離散格式為邊界二階中心差分格式。將定常計算的結果作為非定常數值計算的初始條件，計算時間步長設置為1.15×10-4s，為保證計算結果收斂，計算20個葉輪旋轉周期。

1.3 數值模擬

采用LES對離心泵各個工況進行非定常數值計算，其外特性曲線如圖2所示。離心泵試驗測試系統如圖3所示，其包括測試泵、進出口管路、動力控制部分及各測試儀器等；動力控制部分包括電動機及變頻控制柜；進口管路上安裝進口蝶閥及真空表，出口管路上安裝有壓力表、電磁流量計及出口蝶閥等，通過調節出口閥門的開度來調整泵的運行工況；根據進出口壓差計算泵揚程，同時通過變頻控制柜測試其功率，結合流量及揚程計算泵的效率。通過對比數值計算與試驗結果可以看出，揚程與效率的大渦模擬計算結果略高于試驗測量結果。在設計工況下，揚程計算誤差為3.7%，且在各個工況中最大誤差不超過5%。由于小流量工況下流動復雜，其數值計算效率的誤差較大，隨著流量的增加，計算誤差越來越小，在設計工況下效率誤差為1.8%?？傮w上數值計算結果與試驗結果近似一致，揚程及效率曲線與試驗吻合良好。

圖2 數值模擬與試驗的外特性對比Fig.2 Comparison of pump performance curves between numerical simulation and experiment

圖3 離心泵測試系統Fig.3 Test system of the centrifugal pump

2 動態模態分解

DMD方法通過提取流動中的模態，從而準確描述流動結構，是對Koopman模態分解方法的近似表達，對于線性與非線性系統均可進行分析。將n個非定常流場快照表示為一個矢量集，即

(1)

式中，第i個時刻的流場快照表示為列向量vi，且任意2個快照之間的時間間隔均為Δt。假定流場vi+1可以通過流場vi的線性映射表示

vi+1=Avi

(2)

式中，A為高維流場的系統矩陣，其將時空流場沿時間維度平移Δt，因此矩陣A能夠反映動態系統的時間演化特性。

隨著采集快照數量的增多，基于Arnoldi算法[27]，最后一個流場快照可以表示為之前所有快照的線性組合

(3)

式中：aT=[a1,a2,…,an-1]；r為殘差向量。根據式(4)可以得到

(4)

或其矩陣形式

(5)

(6)

通過對伴隨矩陣C的特征分解，得到其特征值，它是系統矩陣A特征值的近似。

(7)

(8)

因此，矩陣A可近似為

(9)

Φj=Uyj

(10)

gj=Re{lg(λj)}/Δt

(11)

fj=Im{lg(λj)}/Δt

(12)

定義模態振幅α為

α=Y-1UHv1

(13)

式中：Y為特征向量；αi為第i個模態的振幅，代表了該模態對初始流場v1的貢獻。對于標準DMD算法，其模態按照該振幅進行排序，本研究中根據模態范數表示的能量對各階模態進行排序。

根據上述公式，任意時刻的流場可以被重構為

(14)

3 結果與分析

3.1 離心泵葉輪內非定常特性分析

對離心泵全流場進行LES數值計算，并對其葉輪內非定常特性進行分析。葉輪旋轉軸為z軸，從進口方向看沿逆時針方向旋轉。沿z軸正向分別取三個不同截面進行流場對比分析，如圖4所示。S代表與葉輪后蓋板的距離比，其中S=0.5表示葉輪中間截面。

圖4 離心泵葉輪內各截面位置Fig.4 Position of each section in centrifugal pump impeller

圖5為設計工況下葉輪不同軸垂面上相對速度分布，可以看出葉輪內不同軸垂面上速度分布存在較大差別。在S=0.1截面上葉輪內相對速度分布基本均勻，葉片工作面存在明顯的低速區，葉片背面流速由進口到出口逐漸增大，且在葉輪出口存在較大的射流區域。隨著截面逐漸靠近前蓋板，葉片背面的速度梯度越來越大，其流動逐漸復雜，開始發生不同程度的流動分離；葉片工作面的低速區域逐漸減小，葉片背面進口附近流速越來越大，且葉輪出口的尾跡區域逐漸增大。在靠近葉輪前蓋板的截面，即S=0.9，葉輪內流動較為復雜，速度梯度較大，葉片背面出現明顯的流動分離結構。

圖5 設計工況下葉輪內不同截面上相對速度分布Fig.5 Relative velocity distribution on different sections in impeller under the designed flow rate condition

圖6為葉輪內S=0.9截面上不同時刻相對速度的分布情況。從圖6中可以看出，設計工況下葉輪內大部分區域流動基本穩定，流線分布光滑且均勻。但在葉片背面存在尺度較小的周期性分離渦結構。在小流量工況下，即0.6Qd(Qd表示設計工況流量)，葉輪流道內流動復雜，存在周期性演變的大尺度旋渦。針對CH1通道進行分析，T/6～3T/6，通道內兩個尺度較小的旋渦融合成一個大尺度旋渦，其充滿整個CH1通道的進口，對其造成一定的堵塞，4T/6時大尺度旋渦逐漸開始分離，且隨著時間的推移逐漸耗散。

圖6 S=0.9截面上不同時刻相對速度分布Fig.6 Relative velocity distribution at different moments on S=0.9 section

對葉輪葉片背面隨葉輪旋轉的相對運動壓力脈動變化進行監測分析，如圖7(a)所示，在葉片背面設置2個監測點B1及B2，對其壓力數據進行快速傅里葉變換(fast Fourier transform, FFT),得到其頻域變化，如圖7(b)、圖7(c)所示?？梢钥闯鋈~輪內主要頻率與葉輪旋轉頻率(24.17 Hz)近似一致，其余峰值基本為葉輪轉頻的倍數，因此葉輪、蝸殼的動靜干涉是造成離心泵葉輪內不穩定流動最主要的因素。

圖7 葉片背面壓力頻譜圖Fig.7 Pressure spectrum on the suction side of blade

盡管對離心泵葉輪內部各監測點數據的頻域分析能夠得到葉輪內主要頻率及其激勵源，但是其對于不同頻率所對應的流場相干結構無法解耦識別，且對于不同頻率及流場結構下的流動物理機制難以進行詳細分析。因此有必要對其內瞬態流場特征進行分解，分析不同頻率流場結構及機理。

3.2 設計工況下速度場DMD分析

根據3.1節分析，設計工況下模型泵葉輪內流動在靠近前蓋板時較為復雜，因此選取S=0.9截面上流場進行DMD分解，進一步對其非定常流動特性進行深入分析。由于旋轉機械內流場具有周期性，因此選取葉輪內單個流道CH1進行分析，快照數據則選取大渦模擬計算穩定后葉輪旋轉一周內180個流場(每2個時間步長選擇一個流場數據)的相對速度分布數據。

圖8為對設計工況下速度場進行DMD分解得到的特征值分布，從圖中可以看出除少數點外基本上所有點均處于單位圓附近。若某個特征值位于單位圓上，則說明該特征值對應的流場結構處于穩定狀態，若特征值位于單位圓內，表示其對應的流場收斂，若特征值位于單位圓外，則其對應的流場逐漸發散，因此葉輪內速度場DMD模態基本為穩定模態或周期性模態。由于特征值為復數，且共軛成對出現，因此圖8中各模態的特征值關于Im(λ)=0上下對稱分布。圖9為各階模態的能量分布，選取前4階能量較高的模態，其對應頻域分別為0，22.20 Hz，73.65 Hz，127.42 Hz。

圖8 1.0Qd速度場DMD特征值分布Fig.8 DMD eigenvalue distribution of velocity field at 1.0Qd

圖9 1.0Qd速度場模態的能量分布Fig.9 Energy distribution of velocity modes at 1.0Qd

圖10為采用DMD模態重構的葉輪旋轉周期內速度場動態演化過程。葉輪內速度場由前4階DMD模態疊加重構，與原始速度場對比可知前4階模態基本可以完整表達出整個流場的流動結構，對局部較大尺度流場結構能夠精確刻畫，最后一個時刻的流場重構平均相對誤差為1.37%。因此，DMD方法能夠較好地實現流場的低維近似，較為精確地反映流場的主要流動特征，通過DMD模態分析方法能夠清楚地分析離心泵葉輪內復雜流場的非定常特性。

圖10 DMD重構流場與原始流場對比Fig.10 Comparison between reconstructed flow field and original flow field

設計工況下的速度場前4階模態中，1階模態為基本模態，其頻率為0，代表了非定常流場的基本穩態結構，在各階模態中所占能量最高，表征了由于葉輪葉片型線及流動沖角引起的流場特征?；灸B的空間分布如圖11所示，其能夠反映葉輪出口的高速射流區與低速尾跡區及葉片背面距離進口1/3位置的流動分離區域等特征，且隨著時間的變化該模態特征基本穩定。

圖11 1.0Qd速度場基本模態Fig.11 Basic mode of velocity field at 1.0Qd

設計工況下速度場2階模態的頻率與葉輪旋轉頻率一致，反映了由于葉輪旋轉，蝸殼與葉輪的相互作用對葉輪內流場結構的影響。當葉片掃掠蝸殼隔舌時葉輪出口背壓發生變化，使得葉輪流道內流動發生脈動，誘導葉片背面發生流動分離。2階模態演化過程如圖12所示，可以看出葉片背面不斷產生正負交替的高能流體團，在其向下游發展的過程中逐漸拉伸并擴散，且在葉輪出口形成尾跡。在葉輪旋轉周期內，葉片背面共產生一對正負交替且向下游發展的流體團，因此設計工況下2階模態表征了葉輪內流動受蝸殼干擾在葉片背面產生流動分離的流場特征。

圖12 1.0Qd速度場2階模態Fig.12 Second order mode of velocity field at 1.0Qd

圖13為速度場3階及4階模態的演化過程，其頻率分別近似為葉輪旋轉頻率的3倍與5倍。由圖可知，3階與4階模態的尺度依次減小且均小于2階模態，同時流道內高能流體團數量相比2階模態也是依次增加。葉片背面不斷交替產生逐漸拉伸、脫落的高、低速流體團，且隨著時間的推移，在向下游發展的過程中逐漸耗散。因此3階～4階模態為動靜干擾模態的高次諧波行為，表征了由于葉輪旋轉及蝸殼干擾誘導產生的葉片背面不穩定渦結構分離、脫落的流場特征。

圖13 1.0Qd速度場3階及4階模態Fig.13 Third and fourth order modes of velocity field at 1.0Qd

為進一步分析葉輪內不同頻率流動結構的穩定性，將前4階模態對應系數(第i時刻的模態系數為(λj)i-1αj)隨時間的變化規律進行對比分析，如圖14所示。由于1階模態為基本模態，其不隨時間發生變化，因此圖中未給出。由圖可知2階模態系數振幅較高且隨著時間的推移基本穩定，進一步說明離心泵葉輪內流體脈動的最主要原因是蝸殼對葉輪的干擾。3階與4階模態盡管初始能量較高于2階模態，但隨著時間的演化，系數的振幅逐漸減小，因此其反映了流道內脫落的渦結構在逐漸的破碎并耗散。

圖14 1.0Qd速度場DMD模態系數Fig.14 DMD mode coefficients for velocity field at 1.0Qd

3.3 小流量工況下速度場DMD分析

對小流量120 m3/h(0.6Qd)、S=0.9截面上相對速度場進行DMD分析，圖15為0.6Qd速度場DMD特征值分布，可以看出特征值均位于單位圓附近或單位圓內，說明流動為周期性流動。根據各階模態的能量大小同樣選取前4階能量較高的模態進行分析，如圖16所示。

圖15 0.6Qd速度場DMD特征值分布Fig.15 DMD eigenvalue distribution of velocity field at 0.6Qd

圖16 0.6Qd速度場模態的能量分布Fig.16 Energy distribution of velocity modes at 0.6Qd

前4階模態流場分布如圖17所示，其中1階模態頻率為0，為基本模態，同設計工況DMD基本模態一樣，表征了由于葉輪流道幾何形狀引起的穩態流場特征。2階模態頻率為19.71 Hz，與葉輪轉頻近似一致，反映了由于葉輪旋轉及蝸殼干擾引起葉片工作面進口流動不穩定，且隨著時間推移向背面擴散，同時在葉片背面發生流動分離，與進口不穩定流體融合引起流道堵塞的失速流場特征。3階與4階模態頻率分別為2階模態的2倍與3倍，流體團尺度依次減小，數量依次增加，反映了流道內不穩定流體團破碎、耗散特征。圖18為0.6Qd速度場DMD模態的系數變化曲線，可以看出各階模態脈動幅值依次減小，且對比圖14中1.0Qd速度場DMD模態系數可知，小流量工況下葉輪內流體脈動幅值更大。

圖17 0.6Qd速度場前4階模態Fig.17 First four modes of velocity field at 0.6Qd

圖18 0.6Qd速度場DMD模態系數Fig.18 DMD mode coefficients for velocity field at 0.6Qd

4 結 論

(1) 離心泵葉輪在靠近前蓋板時內部流動復雜，葉片背面易發生流動分離。對葉輪內相對速度場進行動態模態分解后其特征值基本分布在單位圓上或圓內，流動具有周期性，且其流場結構可以分解為不同能量及頻率的流動特征，包括基本模態特征、動靜干擾模態及其高次諧波行為流場特征。

(2) DMD方法能夠剝離出反映主要流場特征的主要模態，其基本模態的頻率為零，代表了非定常流場的基本穩態結構，反映出由葉輪流道幾何形狀引起的穩態流場特征。2階～4階動態模態表征了由于葉輪旋轉及動靜干擾在葉輪流道內產生的流動分離、旋渦脫落及耗散流場特征。

(3) 采用主要模態疊加對葉輪旋轉周期內速度場演化過程進行重構，前4階模態基本可以精確地進行流場預測，重構速度場能夠保留原始速度場中重要的動態信息，實現其低維近似，且最后一個時刻的平均相對預測誤差為1.37%。因此，通過DMD方法，可以提取出若干能夠反映葉輪流道內流動特征的主要流場結構，為離心泵內流場穩定性的提升提供一定的理論支撐。