徐勝利，張輝，金婷婷，楊成立，周演

(1. 大連理工大學能源與動力學院，遼寧 大連 116024； 2. 上海電氣燃氣輪機有限公司，上海 200240； 3. 大連理工大學工程力學系，遼寧 大連 116024)

蝸殼是離心泵的重要過流部件之一，其水力性能對離心泵的整體水力性影響顯著.在離心泵輸送流體過程中，蝸殼中的水力損失達到離心泵整體水力損失的一半以上[1]，蝸殼水力性能優(yōu)化設計對提高離心泵的水力效率具有重要意義.傳統(tǒng)的蝸殼一維設計方法簡單方便，目前發(fā)展已經(jīng)十分成熟，主要包括平均速度法、等環(huán)量法、阿基米德螺旋線方程法、結構方形法、不等距方形法等[2].這些一維設計方法需要基于一定的假設條件，如恒定角動量假設、恒定平均速度假設等[3].蝸殼由于受非軸對稱幾何形狀的影響，其內部的瞬態(tài)流場呈非軸對稱分布[4]，特別是受蝸殼隔舌的影響，使得蝸殼入口處的流動參數(shù)的實際分布與理想假設有差別[5].同時，在一些特殊設計條件下，如外部空間復雜、安裝尺寸受限時，傳統(tǒng)蝸殼設計方法難以滿足對空間位置和尺寸限制的要求.因此，發(fā)展一種不依賴于已有設計經(jīng)驗，能夠在外部空間受限條件下進行自動蝸殼優(yōu)化設計的方法顯得十分重要.

近年來，流體拓撲優(yōu)化方法發(fā)展迅速.拓撲優(yōu)化密度法通過多孔介質模型建立材料密度與流動阻力的關系，實現(xiàn)設計域中固體和流體之間的轉化[6].起初，流體拓撲優(yōu)化主要應用在低雷諾數(shù)流動問題的流道設計中[7-8].隨著凍結湍流方法的提出[9]，湍流黏度的變化可以忽略，流體拓撲優(yōu)化也能夠應用到湍流問題的優(yōu)化設計中[10-11].由于拓撲優(yōu)化設計變量多，通常采用梯度優(yōu)化算法求解優(yōu)化問題，伴隨方法是性能函數(shù)靈敏度計算的主要方法.伴隨方法通過原始控制方程和性能函數(shù)導出伴隨方程和伴隨邊界條件，通過求解伴隨方程即可快速得到靈敏度的值.應用連續(xù)伴隨法，已經(jīng)實現(xiàn)了二維、三維管道流動的拓撲優(yōu)化[10].近年來，流體拓撲優(yōu)化方法被應用到離心泵葉輪的流道拓撲優(yōu)化設計中[12-13].離心泵試驗測試也驗證了拓撲優(yōu)化方法能夠提升泵的整體性能[14].

文中采用流體拓撲優(yōu)化方法進行離心泵的蝸殼流道優(yōu)化設計，建立適用于外部空間受限條件下蝸殼設計的拓撲優(yōu)化方法，基于開源軟件OpenFOAM建立相應的優(yōu)化問題分析和求解程序.針對3種不同的蝸殼出口位置，進行蝸殼流道優(yōu)化設計，獲得合理的蝸殼形狀.

1 優(yōu)化問題

1.1 優(yōu)化列式

以離心泵的蝸殼為研究對象，構建了合適的設計域，并假設設計域為多孔介質模型，對其進行優(yōu)化設計，優(yōu)化問題可描述為

find:ρ=(ρ1,ρ2,ρ3,…,ρn),

s.t.:R[v,p,α(ρ)]=0,

(1)

0≤ρi≤1, (1≤i≤n),

式中：ρ為單元密度；下標n為設計域內單元總數(shù)目；J為目標函數(shù)，即設計域進出口的總壓差與流體密度的比值，m5/s3；p為靜壓與流體密度的比值，m2/s2；v為流體速度；S為設計域的邊界；n為邊界面的單位法向量；R為流動控制方程；α為多孔介質孔隙度的倒數(shù)；V為單元體積；m為流體域所占設計域總體積的比例上限.

當單元密度ρ=1時，α=αmax，對應為固體單元；ρ=0時，α=αmin，對應為流體單元.對于ρ在0～1之間的中間密度值，通過基于Heaviside函數(shù)的非線性密度過濾方法，使中間密度值趨于0或1.流體在多孔介質中流動時受到的阻力被稱為Darcy阻力，其大小與孔隙度的倒數(shù)α以及流體流動速度有關，可表示為f=αv.將Darcy阻力項加入不可壓縮、穩(wěn)態(tài)的Navier-Stokes方程中，即可得到文中所使用的控制方程

(2)

(3)

式中：ν為運動黏度，對于湍流，ν可認為是流體的運動黏度與湍流黏度之和.

1.2 密度過濾

在拓撲優(yōu)化中，優(yōu)化結果可能會出現(xiàn)棋盤分布現(xiàn)象，即固體單元ρ=1和流體單元ρ=0交替排列，類似于棋盤狀圖案，不能得到可用的優(yōu)化設計方案.通過非線性密度過濾來控制設計域內的密度分布，能夠消除優(yōu)化過程中可能出現(xiàn)的棋盤分布現(xiàn)象.首先對單元某一特定范圍內的全部單元密度進行加權平均，再以密度平均值作為范圍中心這一單元的線性過濾密度.線性密度過濾如圖1所示，圖中C為線性密度過濾范圍圓心，Rg為過濾半徑.

圖1 線性密度過濾示意圖

(4)

式中：Xi為過濾范圍內單元i中心點的坐標；XC為單元C中心點的坐標；w(Xi)為過濾范圍內所有單元密度的權重，其與Xi的對應函數(shù)關系為

(5)

(6)

由于式(6)為分段函數(shù)，故需要對其進行連續(xù)化處理，如式(7)所示，以此保證靈敏度的連續(xù)性.

(7)

式中：β為控制函數(shù)凹凸程度的參數(shù).圖2為β取不同值的非線性函數(shù)曲線.

圖2 連續(xù)化Heaviside函數(shù)

(8)

式中：αmin為多孔介質孔隙度倒數(shù)α的下限(文中為0)；αmax為α的上限，理論上αmax應趨近于無窮大，在實際計算中選取一個較大的數(shù)來近似.

1.3 伴隨方程

采用梯度類優(yōu)化算法求解優(yōu)化問題，需要計算目標函數(shù)對所有設計變量的靈敏度

(9)

(10)

式中：Ω為整個設計域.拉格朗日乘子可表示為

(u,q)=(u1,u2,u3,q).

(11)

對拉格朗日函數(shù)關于設計變量求全導，可得

(12)

顯然

(13)

目標函數(shù)J與設計變量ρ無關，因此式(13)中等號右端第1項為0.

目標函數(shù)J可以表示為設計域內部積分與邊界積分之和，即

(14)

將方程(2)，(3)代入到式(13)中，結合分部積分公式，整理后可得

(15)

(16)

(17)

在求解控制方程時，采用的入口邊界條件為壓力梯度為0，指定入口速度；出口邊界條件為壓力為0，速度梯度為0；壁面處采用無滑移邊界條件.由此得到伴隨邊界條件.

入口及壁面處

(18)

出口處

(19)

在滿足伴隨方程和伴隨邊界條件之后，式(15)剩下的部分即為靈敏度表達式

(20)

以上的推導過程采用了凍結湍流的方法，忽略了湍流黏度的變化，因此靈敏度求解結果對于層流來說是準確的，而對于湍流流動，得到的靈敏度是近似的.蝸殼優(yōu)化算例結果表明，凍結湍流方法得到的靈敏度具有足夠的精度，能夠滿足蝸殼優(yōu)化的要求.

1.4 優(yōu)化流程

通過改變材料的分布來改變不同單元中流體所受的阻力，流體的流速及壓力分布也隨之有所改變，進而可以影響目標函數(shù)的值.由此，對設計域進行優(yōu)化就可以獲得在一定約束條件下的最優(yōu)材料分布，即蝸殼的形狀.蝸殼的優(yōu)化流程為

1) 對初始設計變量場進行密度過濾;

2) 通過OpenFOAM求解控制方程和伴隨方程，獲得速度、伴隨速度等數(shù)據(jù);

3) 計算目標函數(shù)對于設計變量的靈敏度;

4) 將目標函數(shù)和靈敏度傳遞給優(yōu)化程序，采用MMA優(yōu)化器求解后獲得新的設計變量值并判斷是否收斂，收斂條件為相鄰2次優(yōu)化計算中目標函數(shù)的相對誤差不超過10-5，或優(yōu)化達到最大迭代次數(shù)200次;

5) 若未達到收斂條件，則將原設計變量場更新為MMA優(yōu)化器求解后獲得的設計變量場，并返回第一步重新進行密度過濾及后續(xù)計算，直到達到收斂條件.

2 優(yōu)化算例介紹

2.1 設計域構型

根據(jù)單級單吸清水離心泵的設計參數(shù)構建合適的蝸殼設計域，泵的質量流量Q為55.56 kg/s，揚程H為18.5 m,轉速n為1 480 r/min.

蝸殼的螺旋部分是蝸殼內產生水力損失的主體部位，因此對螺旋部分單獨進行拓撲優(yōu)化設計.采用的二維蝸殼設計域構型如圖3所示.蝸殼的基圓直徑為0.27 m，底邊長及總高度為0.47 m，蝸殼出口位于入口圓的正上方(蝸殼的第8斷面處)，出口總高度為0.1 m，出口中心高度h0為0.05 m.約束條件中流體域的體分比上限m設置為0.45.

圖3 蝸殼設計域構型

另外，在某些外部空間受限條件下，蝸殼的出口位置可能需要一定的升高，因此構建了圖4中的2種蝸殼設計域進行優(yōu)化設計，其出口中心高度h分別設置為0.15 m和0.25 m，體分比上限m分別設置為0.40和0.35，其他幾何參數(shù)不變，以此驗證該蝸殼拓撲優(yōu)化設計方法的靈活性和通用性.

圖4 出口中心高度升高的蝸殼設計域構型

2.2 流體流動狀態(tài)

采用雷諾平均方程(RANS)下的k-ε湍流模型來描述流體的流動狀態(tài).k和ε分別代表湍流動能場和湍流耗散率場，其計算公式如下

(21)

(22)

式中:v′1，v′2，v′3分別為流體速度在x，y，z方向的脈動速度；Cμ為k-ε模型的標準參數(shù)，其值為0.09；l為湍流尺度量.

蝸殼入口速度的軸面分量為1.2 m/s，圓周分量為21 m/s，計算得到k和ε的初始值分別為1.098 m2/s2和18.905 m2/s3.在壁面上，k和ε分別采用kqRwallFunction和epsilonWallFunction邊界條件.

3 優(yōu)化結果

3.1 靈敏度分析結果

采用差分法驗證靈敏度.隨機選取設計域的3個單元分別進行靈敏度的驗證，驗證結果如表1所示.由于求解靈敏度時采用了凍結湍流假設，因此湍流情況下的靈敏度求解結果與差分法結果有一定偏差.

表1 靈敏度驗證結果

3.2 優(yōu)化結果

對蝸殼算例進行拓撲優(yōu)化計算，最終得到蝸殼出口中心高度h分別為h0，3h0和5h0的設計域的優(yōu)化結果如圖5所示.

圖5 不同出口中心高度的蝸殼優(yōu)化結果

3.3 蝸殼傳統(tǒng)設計結果

在同樣的設計參數(shù)下，采用傳統(tǒng)的一元理論蝸殼設計方法[15]進行蝸殼設計.蝸殼的基圓直徑與圖4中蝸殼優(yōu)化設計域的參數(shù)相同，為0.27 m.

根據(jù)離心泵設計參數(shù)確定泵的比轉數(shù)ns

(23)

式中：Qv為泵的體積流量.

計算得泵的比轉數(shù)ns=143，根據(jù)比轉數(shù)與隔舌安放角的關系[15]確定蝸殼隔舌安放角φ0=26.9°，計算得蝸殼的隔舌螺旋角α0=3.27°.

蝸殼斷面形狀采用梯形，根據(jù)速度系數(shù)法確定蝸殼斷面1—8的面積.

蝸殼擴散管的進口近似為蝸殼第8斷面，根據(jù)計算得到擴散管排出口徑Dd=0.13m.根據(jù)擴散管進口當量直徑和排出口徑確定擴散管高度L=0.23 m，并計算得到擴散角θ=7.8°.

根據(jù)以上確定的蝸殼各部分結構參數(shù)，得到二維蝸殼構型如圖6所示.

圖6 傳統(tǒng)設計方法下的蝸殼構型

將傳統(tǒng)設計方法下的蝸殼第8斷面以后的擴散管部分去掉.根據(jù)蝸殼拓撲優(yōu)化結果進行精細建模，兩者采用與蝸殼優(yōu)化設計相同的邊界條件、湍流模型等進行流場及目標函數(shù)值計算.2種蝸殼的構型及速度場對比如下.

圖7 2種設計方法下的蝸殼構型對比

由以上對比可以看出，2種設計方法下的蝸殼構型具有相似性，且通過拓撲優(yōu)化得到的蝸殼流動更均勻，因此目標函數(shù)值與通過傳統(tǒng)設計方法得到的蝸殼相比較小，表明在相同的設計條件下，蝸殼拓撲優(yōu)化設計方法可以替代傳統(tǒng)的蝸殼設計方法，獲得較為合理的蝸殼構型.對于出口中心高度分別為3h0和5h0的蝸殼設計域，傳統(tǒng)的一元理論蝸殼設計方法沒有明確指出蝸殼斷面之間應如何過渡，無法進行此類蝸殼設計；而拓撲優(yōu)化方法能夠很好地適應設計域構型的改變，快速得到符合設計要求的蝸殼構型.

圖8 2種設計方法下的蝸殼速度場對比

4 結 論

1) 將基于連續(xù)伴隨的拓撲優(yōu)化方法應用到了蝸殼的設計中，以總壓損失為目標函數(shù)，推導了所對應的伴隨方程和伴隨邊界條件，在OpenFOAM中實現(xiàn)了用連續(xù)伴隨法求解目標函數(shù)對于設計變量的靈敏度，在一定的體積約束下，得到了二維蝸殼的優(yōu)化構型.從結構來看，用基于連續(xù)伴隨法的流體拓撲優(yōu)化方法得出的蝸殼結構與用傳統(tǒng)設計方法設計出的蝸殼結構十分接近；且拓撲優(yōu)化方法下的蝸殼流場內的總壓損失更小，由此可證明此蝸殼設計方法的可行性；同時也對不同出口中心高度的蝸殼進行了優(yōu)化，最終快速得到了對應的二維蝸殼構型，驗證了此方法的靈活性與通用性.

2) 目前，雖然蝸殼的傳統(tǒng)設計方法已經(jīng)十分成熟，但是對于一些在外部空間受限條件下工作的蝸殼，傳統(tǒng)設計方法設計經(jīng)驗不足.而蝸殼的流體拓撲優(yōu)化方法由于其不依賴于設計經(jīng)驗，且能夠快速得到指定設計條件下的優(yōu)化設計構型，因此在未來的工作中，可以被應用于外部安放空間受限以及其他復雜設計條件下的蝸殼設計中.