流固耦合的多層攪拌反應器振動特性研究

胡效東， 王 灝， 王世飛， 張景坡， 曾慶良

(山東科技大學 機械電子工程學院，山東 青島 266590)

流固耦合的多層攪拌反應器振動特性研究

胡效東， 王 灝， 王世飛， 張景坡， 曾慶良

(山東科技大學 機械電子工程學院，山東 青島 266590)

針對某型號攪拌反應器攪拌主軸因振動斷裂的問題，基于流固耦合理論，采用CFD和FEM相結合的單向耦合穩(wěn)態(tài)的數(shù)值模擬以及經(jīng)驗公式計算方法。在考慮流場作用的條件下分析多種方法得出攪拌器的臨界轉速結果。分析主軸斷裂的原因及影響攪拌器的固有頻率的流體因素，并對主軸結構進行改進。得出結論：主軸驅動電機轉速大于主軸的最低設計轉速是主軸產(chǎn)生振動斷裂的主要原因；攪拌介質的流體力和阻尼作用對攪拌器的固有頻率均有影響；可以將主軸適當加粗以提高其臨界轉速，從而遠離共振范圍并有效解決振動問題；采用流固耦合的方法能夠考慮攪拌器轉動、流體載荷和攪拌槳葉的影響。其計算結果比經(jīng)驗公式更接近實際情況。

攪拌反應器；流固耦合；臨界轉速；CFD

攪拌反應器是化工、石油、食品、醫(yī)藥等行業(yè)的關鍵設備之一。在攪拌、混合、分離等多種工藝生產(chǎn)中具有重要作用。其功能的實現(xiàn)主要依靠于攪拌轉子的轉動。所以對攪拌反應器的振動特性和穩(wěn)定性具有較高的技術要求，對實際中的攪拌反應器的振動特性進行研究具有重要意義。

模態(tài)分析則是進行一切較為復雜的動力特性分析的前提和基礎。為了得到系統(tǒng)較為準確的動力特性分析結果，必須首先得到符合實際情況的系統(tǒng)的模態(tài)分析結果。轉子模態(tài)問題一直備受國內(nèi)外學者的關注。徐斌[1]基于流固耦合理論，采用 SIMPLEC 算法分析斜流泵葉輪的強度，確定了葉輪最大變形位置是在葉片外緣處。高紅斌等[2]分析單級懸臂泵軸-電機系統(tǒng)的動態(tài)響應，建立的集中質量力學模型和扭轉振動方程，可以方便、快捷地求出 IS 型單級單吸式離心泵主軸的固有頻率、臨界轉速和主振型。劉厚林等[3]采用雙向和單向的流固耦合方法分析了導葉式離心泵的強度，以單向流固耦合技術分析了離心泵的強度。荊豐梅等[4]在潮汐水輪機的研究中采用實驗法驗證了單向順序耦合的可靠性。黃文俊等[5]采用有限元法計算風機轉子的臨界轉速。潘旭等[6]對高速轉子槳葉結構的軸流泵的葉片進行振動分析。

本文針對某化工企業(yè)的攪拌反應器主軸因振動過大而斷裂問題。基于流固耦合理論，采用 CFD 方法模擬攪拌反應器內(nèi)部的流體特性，對攪拌反應器進行應力分析，并計算流固耦合作用下攪拌反應器的多階固有頻率，并對比分析攪拌介質的流體力和阻尼作用計算臨界轉速。結合相關的經(jīng)驗公式計算結果，探討主軸斷裂原因及改進措施。

1 動力學方程及分析

對于攪拌反應器這樣一個振動系統(tǒng)，可以假設其系統(tǒng)總共具有n個自由度，則可以得出其系統(tǒng)的振動方程為

(1)

(2)

在考慮水對轉子的影響時，考慮流體波動方程并離散化見式(3)。

(3)

結合式(1)可得式(4)。

(4)

式中：p為流體聲學壓力，M為質量矩陣，C為阻尼矩陣，K為剛度矩陣，F(xiàn)為應力矩陣，F(xiàn)f為流體附加激振力矩陣，Kr為離心應力剛度矩陣，Mf為流體等效質量矩陣，ρf為流體密度，R為流固耦合矩陣，Kr為離心應力剛度矩陣，Kf為流體附加剛度矩陣[7-8]。

在假設系統(tǒng)處于自由振動時，系統(tǒng)運動方程式(2)的解如式(5)

Xi=A(i)sin(ωnit+φi)

(5)

對于振動系統(tǒng)，式(5)中振幅A(i)必有非零解，則必有式(6)。

(6)

2 攪拌反應器的模型建立

2.1 攪拌器模型的建立

攪拌反應器包括傳動裝置、攪拌釜體、攪拌主軸攪拌槳葉。轉子模型可參考為單軸承約束的懸臂梁模型，如圖1所示。攪拌軸直徑D=120 mm，上層四寬葉旋槳攪拌器直徑DJ3= 850 mm；葉片寬度h3=340 mm，中層離心式攪拌器直徑DJ2=800 mm，下層四寬葉旋槳攪拌器直徑DJ1= 800 mm；葉片寬度h1=320 mm。上層、下層槳葉厚度均為10 mm，中層槳葉厚度為5 mm。攪拌軸總長L=3 140 mm，反應釜內(nèi)徑D2=2 600 mm。反應釜直邊s=3 185 mm，反應釜上下封頭均為標準橢圓封頭。

圖1 反應釜整體攪拌結構圖Fig.1 The structure of the stirred-tank reactor

攪拌軸與攪拌器的材料為不銹鋼S30408，密度ρs=7 930 kg/m3，彈性模量E=200 Gpa，泊松比v=0.285。攪拌物料密度ρ=1 000 kg/m3。工況轉速為n=220 r/min。因只研究攪拌反應器在工況下的動力特性，主要考慮流固耦合作用，故簡化攪拌系統(tǒng)模型，流體區(qū)域建模如圖2所示。

圖2 攪拌反應器示意圖Fig.2 The mixing model of the rotor-mediumagitator reaction

2.2 網(wǎng)格劃分

攪拌反應器流體區(qū)域的網(wǎng)格劃分如圖3。整個流體區(qū)域體積較大，采取分塊劃分網(wǎng)格方法。流體區(qū)域網(wǎng)格總單元數(shù)為760 244，總節(jié)點數(shù)為647 719，如圖3所示。固體結構網(wǎng)格單元總數(shù)為28 119，節(jié)點數(shù)為92 478，網(wǎng)格劃分情況如圖4。

圖3 流體域網(wǎng)格劃分示意圖 Fig.3 The grid diagram of fluid

圖4 固體域網(wǎng)格劃分示意圖Fig.4 The grid diagram of solid

3 基于流固耦合攪拌反應器應力分析

3.1 攪拌反應器流場模擬計算

由文獻[8-10]可知，攪拌轉子在實際工況中，所受載荷是處于非定常的。因此受到載荷波動的影響，轉子的固有頻率值會圍繞一恒定值波動變化，但是工況載荷數(shù)值的變化范圍較小，且其波動周期與轉子的轉動周期有關。因此將攪拌轉子視為處于穩(wěn)態(tài)流場的作用下，采用單向耦合的方法。選擇標準k-ε湍流模型，將動區(qū)域的速度與軸面速度統(tǒng)一設置為工況轉速，模型設定為穩(wěn)態(tài)求解。最后進行迭代運算，導出流場對攪拌反應器的應力值，其流場模擬情況如圖5。

3.2 攪拌反應器靜應力分析

將流場分析結果導入靜應力分析模塊，對攪拌反應器的模型施加約束和載荷。載荷包括攪拌反應器的自身重力，攪拌反應器在工況旋轉速度下產(chǎn)生的旋轉預應力，流場在攪拌過程中產(chǎn)生的流體力。約束包括攪拌軸末端的圓柱約束。因軸承和軸封的非線性影響不適合于后期線性模態(tài)分析方法，而且模型中的攪拌反應器屬于單支撐懸臂梁模型，與傳統(tǒng)的轉子動力學的雙軸承模型有較大差別。因此簡化為圓柱約束，忽略軸承的阻尼和剛度對系統(tǒng)的影響。對流場模擬結果進行分析并導入至攪拌反應器的對應面上如圖6。

圖5 流場流速矢量圖Fig.5 The velocity vector of fluid

圖6 耦合面壓力Fig.6 The stress of coupling surface

最后將設置求解總應變和總應力，對攪拌反應器進行受力分析。等效應力云圖如圖7所示，應變量如圖8所示。由圖7和圖8可知，攪拌槳葉在槳葉根部和槳葉末端受力較大。原因是在攪拌過程中，槳葉作為直接的攪拌裝置帶動流體運動，故流體對槳葉的反作用力較大。由于較薄的攪拌槳葉變形積累，在槳葉末端變形較大。并且中層攪拌結構類似離心泵結構，導致流體受力出現(xiàn)明顯的質心偏移，與文獻[10-11]相符。

槳葉與攪拌軸的固定方式為圓柱形約束，這使攪拌軸所受到的橫向力主要來源于槳葉的受力變形。攪拌軸的應力云圖表明，整個攪拌軸的受力變化較為均勻，無明顯應力突變或應力集中，應力最大值不足以導致攪拌軸的斷裂。所以可以排除因應力過大導致攪拌軸失效的原因。因此，考慮攪拌軸的斷裂與攪拌器的振動特性有關，可能會由于工況轉速過于接近其共振頻率產(chǎn)生過大的振動。

圖7 等效應力云圖 Fig.7 Equivalent stress nephogram

4 攪拌反應器的臨界轉速計算

4.1 有限元法的臨界轉速計算

4.1.1 自由狀態(tài)下的臨界轉速計算

首先對轉子做自由振動狀態(tài)的模態(tài)分析，僅對轉子模型施加相應的約束。采用 ANSYS 模態(tài)分析中的線性 Block Lanczos 算法，攪拌轉子在不考慮預應力和介質阻尼(自由振動)前10階的固有頻率如表1。由模擬結果可知，僅前5階段振型如下：第1階的振型為沿x軸的擺動；第2階的振型為沿z軸的擺動；第3階的振型為沿x軸的扭動；第4階的振型為沿z軸的扭動；第5階的振型為沿y軸的扭轉。葉片葉輪的變形量與距離軸線的遠近有關，距離軸線越遠，變形量越大，呈對稱分布。

圖8 應變量圖Fig.8 Deformation

模態(tài)12345頻率/Hz5.3135.31635.4735.4837.22模態(tài)678910頻率/Hz62.7562.8362.9263.0772.18

4.1.2 考慮預應力狀態(tài)下的臨界轉速計算

攪拌轉子的固有頻率是轉子自身的屬性，主要由轉子的質量矩陣和結構剛度矩陣決定。而剛度矩陣會受應力作用而變化。因此以流場模擬的應力計算結果作為模態(tài)分析的邊界條件。并且施加相應的工況載荷如工況轉速下產(chǎn)生的自身旋轉預應力，自身的重力。求解計算后可以得到無阻尼有預應力條件下轉子的前10階固有頻率和振型，如表2。

因模型為低速轉子，受電機和工藝要求的影響，高階固有頻率所轉化的臨界轉速對實際工況的研究意義不大。因此根據(jù)實際情況對計算結果進行篩選，計算所得的一階臨界轉速(可由相應的一階固有頻率轉化)具有參考價值。由表2可知f2=4.55 Hz。只考慮預應力的轉子模態(tài)分析表明：除第3，4，5階以外，其他階固有頻率均有下降。且隨著階數(shù)的升高，預應力降低效果越不明顯，從第1階的14.3%到第10階的3%。因此在計算該低速轉子模型中的一階臨界轉速時，預應力對臨近轉速的計算影響較大。

表2 流固耦合作用下攪拌反應器前十階固有頻率

4.1.3 考慮攪拌介質阻尼作用下的臨界轉速計算

考慮到實際工況中攪拌轉子處于流體介質的影響，其流固耦合作用較復雜。流體介質對轉子的影響也體現(xiàn)在對轉子的阻尼作用。因此結合實際中的攪拌介質范圍，默認介質為不可壓縮無粘度的理想介質。建立攪拌介質-攪拌轉子模型。并通過插入命令的方式實現(xiàn)對攪拌介質邊界條件的施加，即所謂的濕模態(tài)計算。為研究數(shù)值模擬中干濕模態(tài)對攪拌轉子的臨界轉速影響，以有無流體介質的阻尼作用作為差別條件進行對比模擬，可得表3。

表3 阻尼作用下攪拌反應器前十階固有頻率

由表可知，攪拌轉子在攪拌介質的阻尼作用下，各階固有頻率均有明顯的下降。由其是高階固有頻率受攪拌介質的影響較大，前10階的固有頻率均有明顯的下降。從第1階的7.7%到第10階的33.28%，攪拌轉子的固有頻率明顯降低。與文獻[12]中的結果相符合。但該攪拌轉子模型是低速攪拌轉子，由于實際工況要求的限制，同樣取有效計算結果即第一階臨界轉速，即第一階固有頻率f3=4.9 Hz。由此可知，在分析該攪拌轉子的低階臨界轉速時，臨界轉速受阻尼影響較低。但在分析高階的臨界轉速時，流體介質的阻尼作用對高階臨界轉速影響較大。

4.2 經(jīng)驗公式法的臨界轉速計算

根據(jù)《攪拌與混合設備設計選用手冊》中攪拌軸設計可得出攪拌軸的臨界點轉速。通常把工作轉速n低于第一階臨界轉速ncr的軸稱為剛性軸，給定的實際工作轉速為220 r/min，低于仿真分析的第一階臨界轉速273 r/min。這里按照剛性軸進行計算分析。公式中僅考慮了攪拌介質的附加質量從而求出攪拌軸的第一臨界轉速ncr31=306.4 r/min，而在不考慮預應力(當流體默認為無黏性不可壓縮的理想流體時，流體作用可適用于附加質量計算[13])的數(shù)值模擬中得出臨界轉速為ncr0=294.6 r/min。對比理論計算和有限元法仿真分析得到的臨界轉速可看出，兩種方法得到的結果相差較小為2.87%。從而在理論上驗證了有限元法的可靠性，但是比較濕模態(tài)與有預應力的模態(tài)計算結果以及經(jīng)驗公式計算結果，三者仍有一定差值。但計算結果最小值(即為nmin≤0.7ncr=191 r/min)仍低于工況轉速(即n=220 r/min)，即工況轉速恒處于共振轉速區(qū)間，則可知是由于攪拌主軸的結構不合理，臨界振動頻率比較接近電機的驅動頻率，振動過大降低主軸的使用壽命，從而過早失效。轉子的低階固有頻率受攪拌過程產(chǎn)生的預應力和攪拌介質的阻尼作用共同影響。由于本模型中轉子屬于低速轉子，臨界轉速多為一階固有頻率的轉化結果。對比計算結果可知，采用流固耦合計算預應力的方法可以滿足計算需求。

5 優(yōu)化改進

由模擬計算和經(jīng)驗理論公式計算結果可知，實際工況轉速處于攪拌反應器的共振轉速區(qū)間內(nèi)，易發(fā)生共振問題。調節(jié)攪拌軸的軸徑能夠改變主軸的臨界轉速。將攪拌反應器的軸徑由120 mm擴大至150 mm。采用數(shù)值模擬的計算方法得到攪拌器的固有頻率如表4所示。

由表4可得f4=6.0 Hz，臨界轉速ncr4=60×f=360 r/min。該主軸相應的臨界轉速為n4≤0.7ncr4=252 r/min。

表4 優(yōu)化后攪拌反應器前十階固有頻率

再由經(jīng)驗理論公式計算優(yōu)化后的攪拌軸的第一臨界轉速為ncr5=379.25 r/min，響應的主軸驅動頻率最高為n5=265.48 r/min。兩者結果相接近。均遠高于220 r/min的工況轉速。則優(yōu)化設計后的攪拌主軸振動較低，運轉工況良好。

6 結 論

(1)采用數(shù)值模擬計算該攪拌反應器的動力特性，攪拌轉子在流場中的受力分布具有對稱分布的特點，最大應力出現(xiàn)在葉片與輪轂的結合處，在設計中應預防出現(xiàn)應力集中的現(xiàn)象。

(2)在通過數(shù)值模擬的計算轉子的臨界轉速過程中，攪拌介質的流體力和阻尼作用對轉子固有頻率的影響明顯。僅考慮預應力的臨近轉速計算為273 r/min，僅考慮攪拌介質的阻尼作用的仿真結果為297.6 r/min，相應的經(jīng)驗公式計算結果為306.4 r/min。相應的，后兩種方法的計算結果基本一致，相差2.87%，且攪拌轉子的預應力作用對于低階固有頻率影響作用較大，而攪拌介質的阻尼作用對高階的固有頻率影響作用較大。

(3)實際攪拌器驅動電機額定轉速為220 r/min，轉速高于攪拌器設計最低轉速191 r/min，最終造成攪拌反應器攪拌軸振動過大，經(jīng)過1年工作后產(chǎn)生斷裂。通過增大攪拌主軸的軸頸，能夠有效提高攪拌系統(tǒng)的臨界轉速，使之與驅動電機轉速匹配，有效減小攪拌系統(tǒng)的振動，提高產(chǎn)品的使用壽命。在計算該反應釜的攪拌反應器的臨界轉速時，數(shù)值模擬方法充分考慮了流場對攪拌反應器的應力作用，及攪拌槳葉的影響，更能有效的得到攪拌器的實際振動特性。

[1] 徐斌. 混流式水輪機轉輪流固耦合的動力特性分析[D]. 西安: 西安理工大學，2006.

參 考 文 獻

[1] 趙振華，郝曉弘. 局部保持鑒別投影及其在人臉識別中的應用[J]. 電子與信息學報, 2013, 35(2):463-466.

ZHAO Zhenhua, HAO Xiaohong. Linear locality preserving and discriminating projection for face recognition[J]. Journal of Electronics & Information Technology, 2013, 35(2):463-466.

[2] 肖婷，湯寶平，秦毅，等. 基于流形學習和最小二乘支持向量機的滾動軸承退化趨勢預測[J]. 振動與沖擊, 2015, 34(9):149-153.

XIAO Ting, TANG Baoping, QIN Yi, et al. Degradation trend prediction of rolling bearing based on manifold learning and least squares support vector machine[J]. Journal of Vibration and Shock, 2015, 34(9), 149-153.

[3] HE X F, NIYOGI P. Locality preserving projections[C]// Neural Information Processing Systems16. Vancouver: MIT Press, 2004:153-160.

[4] ZHANG Zhenhua, ZHU Xinzhong, ZHAO Jianmin. Image retrieval based on PCA-LPP[C]// 2011 10th International Symposium on Distributed Computing and Applications to Business, Engineering and Science. Wuxi, China, 2011: 230-233.

[5] 王健，馮健，韓志艷. 基于流形學習的局部保持PCA算法在故障檢測中的應用[J]. 控制與決策, 2013, 28(5):683-687.

WANG Jian, FENG Jian, HAN Zhiyan. Locally preserving PCA method based on manifold learning and its application in fault detection[J]. Control and Decision, 2013, 28(5):683-687.

[6] 張曉濤，唐力偉，王平，等. 基于多尺度正交PCA-LPP流形學習算法的故障特征增強方法[J]. 振動與沖擊, 2015, 34(13):66-70.

ZHANG Xiaotao, TANG Liwei, WANG Ping, et al. Fault feature enhancement method based on multiscale orthogonal PCA-LPP manifold learning algorithm[J]. Journal of Vibration and Shock, 2015, 34(13):66-70.

[7] 楊慶，陳桂明，童興民，等. 增量式局部切空間排列算法在滾動軸承故障診斷中的應用[J]. 機械工程學報, 2012, 48(5): 81-86.

YANG Qing, CHEN Guiming, TONG Xingmin, et al. Application of incremental local tangent space alignment algorithm to rolling bearings fault diagnosis[J]. Journal of Mechanical Eigineering, 2012, 48(5):81-86.

[8] JIA Peng, YIN Junsong, HUANG Xinsheng, et al. Incremental Laplacian eigenmaps by preserving adjacent information between data points[J]. Pattern Recognition Letters, 2009, 30(8):1457-1463.

[9] CHIN T J, SUTER D. Out-of-sample extrapolation of learned manifolds[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2008, 30(9):1547-1556.

[10] 談超，關佶紅，周水庚. 基于等角映射的多樣本增量流形學習算法[J]. 模式識別與人工智能, 2014, 27(2):127-133.

TAN Chao, GUAN Jihong, ZHOU Shuigeng. Multi-sample incremental manifold learning algorithm based on isogonal mapping[J].PR & AI, 2014, 27(2):127-133.

[11] 孫即祥. 現(xiàn)代模式識別[M]. 長沙：國防科技大學出版社，2002.

[12] YAN Shuicheng, XU Dong, ZHANG Benyu, et al. Graph embedding and extensions—A general framework for dimensionality reduction[J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2007, 29(1):40-51.

Analysis on a multiple-impeller considering the fluid-structure interaction and its vibration characteristics

HU Xiaodong, WANG Hao, WANG Shifei, ZHANG Jingpo, ZENG Qingliang

(College of Mechanical and Electronic Engineering, Shandong University of Science and Technology, Qingdao 266590, China)

Aiming at analysing the damage mechanism of the stirring spindle of a certain type of mixing reactor due to vibration fault, based on the theory of fluid solid coupling, a numerical simulation, by using the combined method of CFD and FEM, on the one-way coupled stationary vibration was carried out and an empirical formula calculation was also conducted for comparison. Considering the effect of the flow field, the critical speed of the mixer was analyzed. The causes of the spindle fracture and the effect of the fluid factors on the natural frequency of the mixer was analyzed. On this basis, the structure of the spindle was improved. It is concluded that the motor speed, driving the spindle, greater than the minimum design spindle speed is the main reason of the vibration generating the main shaft fracture. The effect of the fluid force and damping on the natural frequency of the mixer is obvious. The enlargement of diameter can increase the critical speed, shift away from the resonance range and solve the vibration problem effectively. The method of analysing the fluid solid coupling can be used to investigate the influences of the rotation of the agitator, the load of the fluid and the mixing blades. The reults are more close to the actual situation than those by the empirical formula.

stirred reactor; fluid solid coupling; critical speed; computational fluid dynamics(CFD)

長江學者和創(chuàng)新團隊發(fā)展計劃(IRT1266 )；國家自然科學基金(51375282)；山東省自然科學基金(ZR2014EEM018)

2016-01-12 修改稿收到日期： 2016-05-31

胡效東 男，博士，副教授，1971年生

E-mail: huxdd@sohu.com

TB533+.1;TF351.5+2

A

10.13465/j.cnki.jvs.2017.14.021