垂直提升管道輸送過(guò)程中的流固耦合效應(yīng)分析＊

周知進(jìn)，盧 浩，王 釗，左明健

（1.湖南科技大學(xué)機(jī)電工程學(xué)院，湖南 湘潭411201；2.湖南科技大學(xué)機(jī)械設(shè)備健康維護(hù)省重點(diǎn)實(shí)驗(yàn)室，湖南 湘潭411201；3.阿爾伯塔大學(xué)機(jī)械工程系，埃德蒙頓，T6G2M7）

流體－結(jié)構(gòu)相互作用問(wèn)題和一般的多物理場(chǎng)問(wèn)題往往過(guò)于復(fù)雜而難以分析解決，所以它們要通過(guò)實(shí)驗(yàn)或數(shù)值模擬方法完成。由于在計(jì)算流體力學(xué)和計(jì)算結(jié)構(gòu)動(dòng)力學(xué)領(lǐng)域的研究工作取得長(zhǎng)足進(jìn)展，這些領(lǐng)域的成果使流體－結(jié)構(gòu)相互作用的數(shù)值模擬得以完成。模擬流體－結(jié)構(gòu)相互作用問(wèn)題有2種主要方法：

單片方法：管流方程和結(jié)構(gòu)位移同時(shí)解決，用一個(gè)單一的求解器一同解決；

分區(qū)方法：分別采用2個(gè)不同的求解方程來(lái)求解流體和結(jié)構(gòu)位移，

單片方法需要這個(gè)物理問(wèn)題特定組合分區(qū)的方法，而保留軟件模塊，因?yàn)橛鞋F(xiàn)成的流動(dòng)求解和結(jié)構(gòu)求解耦合的代碼。此外，分區(qū)的辦法，由于流動(dòng)方程和結(jié)構(gòu)方程的不同，可能更利于開(kāi)發(fā)專門的流體方程和結(jié)構(gòu)方程的解決方案。另一方面，在分區(qū)模擬中需要穩(wěn)定和精確的耦合算法。

同時(shí)Newton－Raphson方法和不動(dòng)點(diǎn)迭代可以用來(lái)解決流固耦合作用（Fluid－Structure Interaction，簡(jiǎn)稱FSI）問(wèn)題。基于Newton－Raphson迭代方法在單片［1－3］和分區(qū)［4－5］方法中被使用。采用 Newton－Raphson的這些方法解決非線性流體方程和結(jié)構(gòu)方程。在Newton－Raphson內(nèi)的迭代線性方程組可以解決系統(tǒng)沒(méi)有雅可比矩陣迭代法先驗(yàn)知識(shí)的系統(tǒng)，而用矢量雅可比產(chǎn)品有限差分進(jìn)行近似。

而Newton－Raphson方法解決在整個(gè)液體和固體域的狀態(tài)流動(dòng)和結(jié)構(gòu)性問(wèn)題，它也有可能重新用于在界面位置未知數(shù)的多自由度系統(tǒng)的FSI裝置問(wèn)題。這個(gè)域分解凝結(jié)成了1個(gè)與界面有關(guān)的子空間FSI誤差［6］問(wèn)題。因此具有未知界面位置的這個(gè)FSI問(wèn)題或者成為了1個(gè)求根或者1個(gè)固定點(diǎn)問(wèn)題。

采用Newton－Raphson迭代的界面牛頓－拉夫森方法解決尋找這個(gè)解的問(wèn)題，例如：雅可比近似線性遞減物理模型［7－8］。在耦合迭代過(guò)程中采用最小二乘模型耦合黑箱子流體界面擬牛頓法雅可比來(lái)逆向逼近求解流體和結(jié)構(gòu)［9］。這種技術(shù)是基于界面擬牛頓最小二乘模型，其中重新表示FSI問(wèn)題作為1個(gè)系統(tǒng)界面的位置和界面上的應(yīng)力分布作為未知數(shù)方程的雅可比矩陣逼近技術(shù)。該系統(tǒng)解決了高斯－賽德?tīng)栴愋偷牧黧w求解和結(jié)構(gòu)求解雅可比矩陣塊擬牛頓迭代逼近最小二乘模型［10］。

固定點(diǎn)問(wèn)題也可能用固定點(diǎn)迭代解決，也被稱為高斯－賽德?tīng)柕?］，這意味著流體問(wèn)題和結(jié)構(gòu)性問(wèn)題都相繼解決，直到變化小于收斂準(zhǔn)則。然而，迭代收斂速度慢，尤其是當(dāng)流體和結(jié)構(gòu)之間的相互作用很強(qiáng)時(shí)，特別是高密度流體／結(jié)構(gòu)比例或不可壓縮的流體的情況［12］。適用基于以前的固定點(diǎn)迭代能被Aitken松馳和最速下降松弛因子穩(wěn)定和加速收斂［10］。

如果流體和結(jié)構(gòu)之間相互作用較弱，僅僅在每1個(gè)時(shí)間步需要1個(gè)固定點(diǎn)迭代。這些所謂的交錯(cuò)或松散耦合方法不強(qiáng)制執(zhí)行1個(gè)時(shí)間步內(nèi)的流體結(jié)構(gòu)界面上的平衡，但它們適用于帶有重型的和相當(dāng)剛性的結(jié)構(gòu)模擬。部分研究分析了用于模擬流體－結(jié)構(gòu)相互作用分割算法的穩(wěn)定性［11－12］。

垂直懸臂提升管道在輸送固液兩相流過(guò)程中，管道（固體）為彈性體，受到內(nèi)外流體作用，流體的作用引起管道壁發(fā)生變形或運(yùn)動(dòng)，管道的變形或運(yùn)動(dòng)又反過(guò)來(lái)改變流場(chǎng)形態(tài)，從而改變流動(dòng)狀態(tài)，流動(dòng)狀態(tài)的改變又會(huì)影響管道的運(yùn)動(dòng)和變形，管道與流體之間的這種相互作用，在不同約束支撐下將產(chǎn)生多種形態(tài)各異的流固耦合現(xiàn)象［3］，即管道與流體之間的耦合作用。以往學(xué)者大多對(duì)水平輸送、固定支撐（多約束）、彈性支承等管道振動(dòng)進(jìn)行研究［4－7］，對(duì)于垂直單端支撐提升管道振動(dòng)研究較少，更沒(méi)有考慮柔性支撐情況下輸送過(guò)程的流固耦合作用，所得結(jié)果精度不高［13－14］。

1 揚(yáng)礦硬管流固耦合數(shù)學(xué)模型

對(duì)揚(yáng)礦硬管工作狀態(tài)進(jìn)行計(jì)算分析的時(shí)候，揚(yáng)礦硬管內(nèi)部固液兩相流體和外部海水必然都是相互影響的。流體作用力施加到管壁（結(jié)構(gòu)）上，結(jié)構(gòu)的變形反過(guò)來(lái)影響流體區(qū)域，這樣就形成了流固耦合過(guò)程。

流固耦合邊界位移協(xié)調(diào)方程［15］：

流固耦合邊界力平衡方程：

式中：df為流體位移；da為結(jié)構(gòu)位移；τf為流體的應(yīng)力；τs為結(jié)構(gòu)的應(yīng)力。

將這些值定義在流固耦合界面上。由位移協(xié)調(diào)方程來(lái)決定流固耦合邊界上流體節(jié)點(diǎn)的位置，另外由力平衡方程決定在流固耦合邊界上的應(yīng)力，根據(jù)公式：

式中：hd為結(jié)構(gòu)節(jié)點(diǎn)的位移；τf為流體的應(yīng)力。流體均布力積分為集中力施加到結(jié)構(gòu)節(jié)點(diǎn)上。

流固耦合中流體和結(jié)構(gòu)模型的時(shí)間積分必須是相容的，在ADINA中流體和結(jié)構(gòu)在流固耦合界面上2個(gè)坐標(biāo)系是相同的，都是用lagragian坐標(biāo)系。流體和結(jié)構(gòu)的方程分別用Gf［f，f′］＝0和Gs［d，d′，d″］＝0來(lái)表示，以f和d分別表示流體和結(jié)構(gòu)的變量。流體的速度和加速度分別表示為：

方程中t＋Δt時(shí)刻的速度和加速度可以用位移未知量來(lái)表示：

將這個(gè)公式應(yīng)用到耦合系統(tǒng)中，可以得出時(shí)間積分格式為：

Gs（t＋Δt）≈Gs［dt＋Δt，dt＋Δta＋ξt，dt＋Δtb＋ηt］＝0 （9）流固耦合有限元方程則表示為：

式中：Xf、Xs分別為流體和結(jié)構(gòu)節(jié)點(diǎn)上的解向量；Ff、Fs分別是Gf和Gs相應(yīng)的有限元方程。

在這個(gè)問(wèn)題中，內(nèi)外部流場(chǎng)作用力引起了結(jié)構(gòu)的變形，同時(shí)結(jié)構(gòu)的位移又影響流場(chǎng)的狀態(tài)，兩者相互作用。在考慮2種作用時(shí)流固耦合分析稱為雙向耦合，流固耦合分析中流固耦合方程必然是非線性的，所以使用迭代方法得到解。ADINA中使用應(yīng)力、位移或者應(yīng)力和位移來(lái)檢查迭代收斂性。

2 揚(yáng)礦管道有限元模型

2.1 管道模型的建立

管道幾何模型采用ADINA—native建模方式，本文根據(jù)不同截面尺寸建立管道幾何模型，按照相關(guān)要求定義和施加邊界條件，一種是對(duì)管道上端按照固接情況施加約束，另一種是對(duì)管道上端按照鉸接情況施加約束，施加載荷類型為Mass Proportional，在Z軸負(fù)方向定義 Magnitude為9.8m／s2，管道單元采用3D solid單元，單元?jiǎng)澐植捎?節(jié)點(diǎn)3D實(shí)體單元，網(wǎng)格密度在線、面、體建模時(shí)依次設(shè)置，并選擇小變形假設(shè)。

2.2 流體模型的建立

將管內(nèi)流體簡(jiǎn)化成圓柱體，采用native建模方式，流體單元采用3DFluid單元，單元?jiǎng)澐植捎?節(jié)點(diǎn)3D實(shí)體單元，并根據(jù)輸送流體體積濃度等情況不同設(shè)置管內(nèi)液體參數(shù)；管外流體區(qū)域大小選擇為管道結(jié)構(gòu)的7倍［15］，單元選擇與劃分和管內(nèi)流體相同。

得到的管道流體流固耦合有限元模型及管道、管外流體、管內(nèi)流體有限元模型見(jiàn)圖1。

圖1 管道、流體、管道流體有限元模型Fig.1 Finite element model of pipe and fluid and pipefluid coupled

3 仿真及結(jié)果分析

3.1 不同輸送礦石體積濃度對(duì)管道固有頻率的影響

揚(yáng)礦管道輸送不同體積濃度礦石顆粒會(huì)對(duì)管道產(chǎn)生不同程度的振動(dòng)，管道固有頻率受礦石顆粒體積濃度影響程度不一樣。

選取同一管道結(jié)構(gòu)，分析不同礦石顆粒體積濃度對(duì)管道固有頻率的影響，選取管段長(zhǎng)度11m，管內(nèi)徑為206mm，壁厚10mm，彈性模量210GPa，泊松比0.3，密度為7800kg／m3。管外流體為水，密度為1 000kg／m3，體積模量為2×109N／m2，流速設(shè)定為0.1m／s。管內(nèi)輸送礦石顆粒密度為2 040kg／m3，體積濃度分別為5%、10%、15%，上端采用固接約束，對(duì)管道進(jìn)行流固耦合模態(tài)分析，計(jì)算3種輸送濃度下管道的前6階固有頻率與相應(yīng)振型，并將管內(nèi)流體—管道二者耦合、管內(nèi)流體—管道—管外流體三者耦合與不考慮流固耦合效應(yīng)進(jìn)行對(duì)比，且定義考慮流固耦合與不考慮流固耦合的管道固有頻率之比為影響系數(shù)［15］。管內(nèi)流體—管道二者耦合固有頻率、管內(nèi)流體—管道—管外流體三者耦合固有頻率仿真結(jié)果如表1、表2所示，相應(yīng)的曲線見(jiàn)圖2～3。

表1 不同輸送濃度的管內(nèi)流體—管道二者流固耦合管道固有頻率Table 1 Natural frequency of fluid－pipe coupled under different transporting volume concentrations

表2 不同輸送濃度的管內(nèi)流體—管道—管外流體三者流固耦合管道固有頻率Table 2 Natural frequency of internal fluid－pipe－outer coupled under different transporting volume concentrations

圖2 不同濃度下管內(nèi)流體－管道耦合頻率Fig.2 Both coupled frequency under different concentration

圖3 不同濃度下管內(nèi)流體－管道－管外流體耦合頻率Fig.3 Three coupled frequency under different concentration

由表1、2可以看出：（1）管道的前6階固有頻率，第1、2階，第3、4階，第5、6階數(shù)值相同，振型也相同，只是振動(dòng)方向不同，即存在強(qiáng)的對(duì)稱性，1、2階振型表現(xiàn)為橫向偏移，3、4階振型表現(xiàn)為單向彎曲，5、6階振型表現(xiàn)為雙向彎曲；（2）管道流固耦合固有頻率隨著輸送濃度的增大而降低，同一階固有頻率的影響系數(shù)也呈遞減趨勢(shì)，即管內(nèi)輸送流體濃度越大，考慮流固耦合的固有頻率較不考慮耦合時(shí)差異越明顯；（3）不同階次固有頻率影響系數(shù)不同，即不能單純的乘以1個(gè)系數(shù)來(lái)處理，在實(shí)際的管道輸送中應(yīng)該考慮流固耦合效應(yīng)對(duì)固有頻率的影響，盡量使管道固有頻率遠(yuǎn)離工作頻率，避免管道共振情況的發(fā)生；（4）考慮管內(nèi)流體—管道—管外流體三者流固耦合時(shí)管道固有頻率較只考慮管內(nèi)流體—管道二者耦合要低，比不考慮流固耦合時(shí)更低，影響系數(shù)也呈下降趨勢(shì)。

從圖2和圖3可以看出，在體積濃度較低情況下，濃度小幅波動(dòng)對(duì)流固耦合管道固有頻率影響較小?？紤]流固耦合與不考慮流固耦合情況隨著模態(tài)階次提高，管道固有頻率差異越來(lái)越大。如輸送濃度5%時(shí)管內(nèi)流體－管道二者耦合時(shí)其固有頻率相差分別為：一階相差0.4Hz、三階相差2.5Hz、五階相差11.15Hz。只考慮管內(nèi)流體－管道耦合的管道固有頻率比考慮管內(nèi)流體－管道－管外流體三者耦合求得的固有頻率高4Hz左右。

3.2 不同管道截面對(duì)管道流固耦合固有頻率的影響

管道的截面尺寸由直徑和壁厚決定，揚(yáng)礦管道根據(jù)設(shè)計(jì)要求具有多種截面，本文用壁厚與管徑之比［14］即管道相對(duì)壁厚d來(lái)討論其對(duì)管道流固耦合固有頻率的影響。下面以管內(nèi)流體為5%體積濃度的礦石顆粒，支撐條件為管道上端固接，計(jì)算各種截面參數(shù)下管道的前6階流固耦合固有頻率。由于其1階與2階，3階與4階，5階與6階固有頻率相同，故選擇1階、3階、5階頻率為對(duì)象討論管道截面尺寸對(duì)流固耦合固有頻率的影響。計(jì)算結(jié)果如表3所示，相應(yīng)曲線如圖4、5、6、7所示。

表3 不同截面參數(shù)下管道流固耦合固有頻率Table 3 Natural frequency of fluid－pipe interaction under different cross－section

從表3和圖2可以看出，流體對(duì)管道流固耦合固有頻率的影響程度隨著相對(duì)厚度d值的增大而減小，因?yàn)閐值越大說(shuō)明管壁相對(duì)越厚，流體對(duì)管道的作用效應(yīng)相對(duì)就弱，則流固耦合效應(yīng)對(duì)固有頻率的影響就小?？紤]管內(nèi)流體—管道—管外流體三者耦合比考慮管內(nèi)流體—管道二者耦合時(shí)各階固有頻率均降低，三者耦合時(shí)相差比二者耦合時(shí)大，說(shuō)明了考慮三者耦合時(shí)計(jì)算結(jié)果更加趨于精確。

4 結(jié)論

本次調(diào)查基于有限元軟件ADINA－FSI模塊的利用，建立了管道流固耦合有限元模型，對(duì)管道－流體之間相互作用進(jìn)行了流固耦合模態(tài)分析，討論了不同輸送體積濃度、不同管道截面尺寸下流固耦合效應(yīng)對(duì)管道固有頻率的影響。結(jié)果表明：

（1）揚(yáng)礦管道流固耦合固有頻率隨著管內(nèi)輸送礦石顆粒濃度的增大而降低，考慮管內(nèi)流體—管道—管外流體三者流固耦合時(shí)的管道固有頻率較只考慮管內(nèi)流體—管道二者耦合和不考慮流固耦合時(shí)要小，影響系數(shù)也相應(yīng)降低，即考慮三者流固耦合時(shí)管道固有頻率結(jié)果更加精確。

圖7 相對(duì)誤差%與管道相對(duì)壁厚b的關(guān)系曲線Fig.7 The curve of relative proportion of FSI natural frequency of pipe with d

（2）流體對(duì)管道流固耦合固有頻率的影響程度隨著相對(duì)厚度d值的增大而減小，在考慮管內(nèi)流體—管道—管外流體三者流固耦合時(shí)管道固有頻率較只考慮管內(nèi)流體—管道二者耦合相對(duì)誤差增大，說(shuō)明了考慮三者流固耦合的重要性。

（3）在實(shí)際的揚(yáng)礦管道輸送中應(yīng)該考慮輸送流體顆粒濃度、不同相對(duì)管道壁厚對(duì)管道固有頻率的影響及管內(nèi)流體—管道—管外流體三者流固耦合作用對(duì)固有頻率的影響，盡量使管道固有頻率遠(yuǎn)離工作頻率，避免共振情況的發(fā)生。

［1］ Heil M.An efficient solver for the fully coupled solution of largedisplacement fluid－structure interaction problems ［J］.Computer Methods in Applied Mechanics and Engineering，2004，193：1－23.

［2］ Bathe K J，Zhang H.Finite element developments for general fluid flows with structural interactions ［J］.International Journal for Numerical Methods in Engineering，2004，60：213－232.

［3］ Hron J，Turek S.A monolithic FEM／multigrid solver for ALE formulation of fluid－structure interaction with application in biomechanics［J］.Lecture Notes in Computational Science and Engineering.Fluid Structure Interaction Modelling，Simulation，Optimisation Springer－Verlag，2006（5）：146－170.

［4］ Degroote J，Bathe K J，Vierendeels J.Performance of a new partitioned procedure versus a monolithic procedure in fluid structure interaction［J］.Computers and Structures，2009，87：793.

［5］ Vierendeels J，Lanoye L，Degroote J，et al.Implicit coupling of partitioned fluid－structure interaction problems with reduced order models［J］.Computers and Structures，2007，85（11－14）：970－976.

［6］ K ttler U，Wall W.Fixed－point fluid－structure interaction solvers with dynamic relaxation［J］.Computational Mechanics，2008，43（1）：61－72.

［7］ Degroote J，Bruggeman P，Haelterman R，et al.Stability of a coupling technique for partitioned solvers in FSI applications［J］.Computers and Structures，2008，86（23－24）：2224－2234.

［8］ Jaiman R，Jiao X，Geubelle P，et al.Conservative load transfer along curved fluid－solid interface with non－matching meshes ［J］.Journal of Computational Physics，2006，218（1）：372－397.

［9］ Vierendeels J，Dumont K，Dick E，et al.Analysis and stabilization of fluid－structure interaction algorithm for rigid－body motion［J］.AIAA Journal，2005，43（12）：2549－2557.

［10］ Vierendeels J，Lanoye L，Degroote J，et al.Implicit coupling of partitioned fluid－structure interaction problems with reduced order models［J］.Computers and Structures，2007，85 （11－14）：970－976.

［11］ K ttler U，Wall W.Fixed－point fluid－structure interaction solvers with dynamic relaxation ［J］.Computational Mechanics，2008，43（1）：61－72.

［12］ Degroote J，Bruggeman P，Haelterman R，et al.Stability of a coupling technique for partitioned solvers in FSI applications［J］.Computers and Structures，2008，86（23－24）：2224－2234.

［13］ 張學(xué)志，黃維平，李華軍.考慮流固耦合時(shí)的海洋平臺(tái)結(jié)構(gòu)非線性動(dòng)力分析 ［J］.中國(guó)海洋大學(xué)學(xué)報(bào)：自然科學(xué)版，2005，35（5）：823－826

［14］ 楊超，范士娟.輸液管道流固耦合振動(dòng)的數(shù)值分析 ［J］.振動(dòng)與沖擊，2009，28（6）：2148－2157.

［15］ 岳戈.ADINA流體與流固耦合功能的高級(jí)應(yīng)用 ［M］.北京：人民交通出版社，2010.