加肋圓柱殼結構的FE-BE算法網格尺度劃分原則研究

黃振衛，周其斗，紀 剛，孟慶昌，段嘉希

(海軍工程大學 艦船工程系，武漢 430033)

加肋圓柱殼結構的FE-BE算法網格尺度劃分原則研究

黃振衛，周其斗，紀 剛，孟慶昌，段嘉希

(海軍工程大學 艦船工程系，武漢 430033)

研究結構有限元耦合流體邊界元算法(FE-BE算法)中結構濕表面的網格尺度劃分原則具有十分重要的意義。以有限長加肋圓柱殼為研究對象，采用波數譜分析方法分離出結構主振型分量(攜帶能量最多的振動模式)。根據分析結果提出FE-BE算法中結構濕表面的網格尺度劃分原則，并采用FE-BE算法驗證了該原則的有效性。

加肋圓柱殼； 結構有限元耦合流體邊界元； 振動與聲輻射； 網格尺度

結構有限元耦合流體邊界元(FE-BE)算法一直被認為是解決水下大型結構低頻振動與聲輻射的有力工具[1-3]，因為采用邊界元方法求解無限域的聲場問題時，邊界元只需要在結構濕表面(結構物與水接觸的表面)上劃分，當結構有限元的劃分完成后，邊界元劃分業已完成。為了節約FE-BE算法的計算成本，一般應選取單元尺度較大的有限元單元劃分結構濕表面，但是這必然帶來計算精度的降低。另一方面，為了保證FE-BE算法的計算精度，結構濕表面的單元尺度應該盡量小以便于描述最短的結構波長，但是這將導致大型復雜結構的邊界元求解規模急劇地增加，使得FE-BE算法的計算時間非常長。因此，對于FE-BE算法，結構濕表面單元尺度的選取與計算精度和計算成本之間有個權衡的問題。目前很少有文獻研究FE-BE算法中結構濕表面的網格尺度劃分原則的問題。

結構的振動可以分解為一系列不同波長的振動分量的疊加，其中攜帶振動能量最大的分量即為結構的主振型分量。一般而言，采用有限元單元劃分結構濕表面時，單元尺度必須要能夠描述結構濕表面上最短的彎曲波長[4]，但是實際振動結構的最短彎曲波長可以無限小，如果采用非常小的單元尺度劃分結構濕表面，將會使得系統的邊界元求解規模急劇地增加。顯然，不能簡單地采用結構最短彎曲波長作為FE-BE算法中濕表面網格劃分的參考標準。由于結構主振型分量的彎曲波長對于流固耦合計算精度有重要的作用，因此，本文重點采用波數譜方法分離有限長圓柱殼的主振型分量，根據分析結果提出FE-BE算法中結構濕表面的網格尺度劃分原則，并采用FE-BE算法驗證了該原則的有效性。本文首次提出了以結構主振型分量的彎曲波長作為FE-BE算法中結構濕表面的網格尺度劃分的參考標準，結論對于FE-BE算法的有限元模型的網格劃分具有十分重要的參考意義。

1 模型及有限元劃分

以有限長加肋圓柱殼為研究對象，結構幾何尺寸及材料和流體屬性如表1所示。殼體內表面上均勻布置23根肋骨，激勵力作用于中間肋骨上，沿著徑向指向殼體外部的方向。采用有限元軟件PATRAN進行建模，結構濕表面采用三角形單元進行有限元劃分，肋骨采用四邊形單元進行有限元劃分，加肋圓柱殼的有限元模型示意圖如圖1所示。72個計算聲場點均布于與柱殼中間肋骨同心，半徑為9 m的圓周上。計算頻率為50 Hz～500 Hz，頻率間隔為5 Hz，模型為自由邊界條件，不考慮自由水面對結構振動與聲輻射的影響。

表1 模型幾何尺寸及材料和流體屬性Tab.1 Geometric size and parameter values for calculation

圖1 加肋圓柱殼有限元模型Fig.1 Finite element model of the ring stiffened cylindrical shell

2 主振型的分離

本節采用編制的圓柱殼結構振動波數譜程序，分離柱殼結構主振型分量，該程序的正確性已在文獻[5]中得到很好的驗證。為了描述有限長圓柱殼結構做剛體運動的波數譜，將柱殼長度延長為原柱殼長度的12倍，假定延長的部分的外殼體法向速度為0[5]，如圖2示，L為圓柱殼的長度，a為圓柱殼的半徑，圖中省略了肋骨示意。

圖2 無限長柱殼有限長模型Fig.2 An infinite length cylindrical shell with finite vibration distribution

采用編制的FE-BE算法程序計算圓柱殼濕表面的法向位移復數幅值w(z,θ)，即

(1)

式中：wR和wI分別為是結構濕表面法向位移的實部和虛部。值得注意的是，本節為了保證FE-BE算法程序計算結構濕表面振動響應的精度，加肋圓柱殼的每檔肋骨間距采用6個單元進行劃分，周向采用48個單元進行劃分，最終得到15 552個結構單元，8 330個節點，14 448個濕表面單元，7 226個濕表面節點。

將法向位移的實部wR沿著柱殼的軸向和周向兩個方向進行傅里葉變換，得到實部所代表的駐波場

(2)

式中，kz為軸向波數，n為周向完整波的數量，且有

(3)

(4)

同理，可以得到虛部所代表的駐波場

(5)

總的法向速度的振動功率可表示為

(6)

式中:Ev(n,kz)為法向速度振動功率譜，表達式為：

Ev(n,kz)=

(7)

根據式(7)，即可繪制出不同的周向振動模式n時，結構法向速度振動功率隨著結構軸向波數的變化曲線，進而分離出加肋圓柱殼的主振型分量以及主振型分量的彎曲波長λb。

圖3為模型的法向速度振動功率譜(n=0,1,…,9)，從圖3(a)中可以看出，當n=2，無因次軸向波數大約為0時，結構波的振動能量最大(主振型分量)，對應的結構軸向波長為∞。分析圓柱殼周向方向上結構波的傳播特性，可以得到n=2時，周向上的結構波長為2πa/n=1.99 m。值得一提的是結構中其它非主振型分量攜帶的能量小于主振型分量所攜帶的能量，例如圖3中振動能量級大3 dB，意味著攜帶的能量將增大一倍左右。

采用相同的方法即可分離出其他頻率下結構主振型分量，計算結果見表2。

(a) 50 Hz

(b) 300 Hz

(c) 500 Hz

表2 不同頻率下結構主振型分量Tab.2 Main vibration mode at different frequency

3 網格尺度對FE-BE算法的影響

由于加肋圓柱殼的振動與聲輻射缺乏理論解，本文將逐步加密結構濕表面上的網格，觀察FE-BE算法的計算結果是否趨于穩定，并與NASTRAN軟件中的FE-IE算法[6]計算結果進行對比。為了保證無限元法的計算精度，在結構外域構造一個長10 m，半徑5 m的圓柱面人工邊界，人工邊界的單元尺度為0.483 m，結構濕表面軸向單元尺度為0.059 m，周向單元尺度為0.083 m，無限元的計算階數為10階。根據第2節波數譜的計算結果，設計工況如表3所示。為了便于研究，所有工況下柱殼圓周方向上的單元尺度選取為Lc=0.083 m，定義Nz=λb/Lz為一個主振型波長內的單元數量，其中，Lz為結構濕表面軸向方向上的單元尺度，λb為結構主振型分量軸向彎曲波波長。

表3 計算工況Tab.3 Computational cases

圖4～圖6為各工況下輻射聲壓級的FE-BE算法計算值對比，從圖中可以看出：

(1) 激勵頻率為50 Hz時，不同結構濕表面網格尺度下FE-BE算法的計算結果均與FE-IE算法的計算結果吻合良好。主要原因是，此時結構主振型的軸向彎曲波波長為無窮大，即使軸向單元尺度選取為一個肋骨間距一個單元(Lz=0.238 m)時，也能準確描述主振型分量的彎曲波波長。但是結構濕表面的單元尺度(即邊界元尺度)不能無限增大，根據文獻[7]，邊界元的單元尺度一般應保證一個聲波長至少有6個單元，因此，結構濕表面的單元尺度最大不能超過1/6倍聲波長；

(2) 激勵頻率為300 Hz時，隨著結構濕表面的單元尺度減小，FE-BE算法的計算結果趨于穩定，且與FE-IE算法的計算結果吻合良好。主要原因是300 Hz時FE-BE算法的結構濕表面網格尺度保證了一個主振型波長內至少有7個有限元單元，足以描述結構的主振型分量的波動特性；

(3) 激勵頻率為500 Hz時，隨著結構濕表面的單元尺度減小，FE-BE算法的計算結果趨于穩定，但與FE-IE算法的計算結果有一定的差別。當結構濕表面的軸向單元尺度為Lz=0.238 m時，一個主振型彎曲波長內只有3.2個單元，無法準確地描述結構主振型分量的波動特性，因此，單元尺度為Lz=0.238 m的FE-BE算法計算結果與濕表面網格加密后FE-BE算法計算結果差別較大。當結構濕表面的軸向單元尺度為Lz=0.119 m時，保證了一個主振型彎曲波長內有6.5個有限元單元，此時FE-BE算法的計算結果與濕表面網格加密后FE-BE算法計算結果吻合良好。此外，500 Hz時網格加密后FE-BE算法計算結果與FE-IE算法計算結果產生一定差別的原因可能是近場效應對FE-IE算法的影響更為明顯，為了改善FE-IE算法的計算精度，人工邊界應該選取在距離結構濕表面更遠的位置或者采用更小的單元尺度劃分人工邊界。

圖4 加肋圓柱殼在50 Hz時的聲壓級Fig.4 Numerical SPL of the ring stiffened cylindrical shell at 50 Hz

圖5 加肋圓柱殼在300 Hz時的聲壓級Fig.5 Numerical SPL of the ring stiffened cylindrical shell at 300 Hz

圖6 加肋圓柱殼在500 Hz時的聲壓級Fig.6 Numerical SPL of the ring stiffened cylindrical shell at 500 Hz

總的來說，對加肋圓柱殼結構而言，當結構主振型分量的彎曲波波長小于對應激勵頻率下的聲波波長時，以結構主振型分量的彎曲波長作為FE-BE算法中結構濕表面的網格尺度劃分標準是可行的，即結構濕表面的網格尺度應保證一個主振型分量的彎曲波長內至少有6個單元。值得注意的是，此時結構濕表面的單元尺度能夠描述結構中攜帶能量最大的振動分量，進而保證了FE-BE算法的計算精度。若將結構濕表面的單元尺度進一步減小或者增大，都有可能帶來對應的計算成本的增加或計算精度的損失。

為了進一步應用FE-BE算法的網格劃分原則，本節將計算50～500 Hz時外域聲場中A(0 m,100 m,0 m)點的輻射聲壓級。根據以上分析，計算頻率在500 Hz以內時，加肋圓柱殼每個肋骨間距采用2個單元進行劃分(Lz=0.119 m)即可保證FE-BE算法的計算精度。同時，FE-IE算法中加肋圓柱殼每個肋骨間距4個有限元單元進行劃分，人工邊界為一個長12 m，半徑為8 m的圓柱面，人工邊界的單元尺度為0.48 m，無限元的計算階次為10階。A點的輻射聲壓級傳遞函數頻響曲線如圖7所示，從圖中可以看出，激勵頻率較低時，FE-BE算法的輻射聲壓級傳遞函數頻響曲線與FE-IE算法的頻響曲線吻合良好；激勵頻率較高時，FE-BE算法的輻射聲壓級傳遞函數頻響曲線與FE-IE算法的頻響曲線變化趨勢一致。

圖7 結構聲學傳遞函數的FE-BE計算值Fig.7 Sound transfer function obtained by the FE-BE method

4 結 論

本文重點采用波數譜方法分離圓柱殼結構的主振型分量，根據分析結果提出FE-BE算法中結構濕表面的網格尺度劃分原則，并采用FE-BE算法驗證了該原則的有效性，計算結果表明：對于文中算例而言，當結構主振型分量的彎曲波長小于對應激勵頻率下的聲波長時，以結構主振型分量的彎曲波長作為結構濕表面的網格尺度劃分標準是可行的，即結構濕表面的網格尺度應保證一個主振型分量的彎曲波長內至少有6個單元。

[1] ZHOU Q, JOSEPH P F.A numerical method for the calculation of dynamic response and acoustic radiation from an underwater structure[J].Journal of Sound and Vibration,2005(283): 853-873.

[2] 紀剛,張緯康,周其斗.靜水壓力作用的水下結構振動及聲輻射[J].中國造船,2006,47(3):37-44.

JI Gang, ZHANG Weikang, ZHOU Qidou.Vibration and radiation from underwater structure considering the effect of static water preload[J].Shipbuilding of China,2006,47(3):37-44.

[3] EVERSTIN G C, HENDERSON F M.Coupled finite element/boundary element approach for fluid-structure interaction[J].J.Acoust.Soc.Am., 1990,87:1936-1947.

[4] 宗福開.波傳播問題中有限元分析的頻散特性及離散化準則[J].爆炸與沖擊,1984,3(4):16-23.

ZONG Fukai.Frequency dispersion characteristic and discretization of the finite element analysis in wave propagation problems[J].Explosion and Shock Waves,1984,3(4):16-23.

[5] 譚路.基于波數譜的細長殼體圓柱殼振動與聲輻射噪聲特性研究[D].武漢：海軍工程大學，2013.

[6] 田麗思,李相輝,馬越峰,等.MSC Nastran 動力分析指南[M].北京：中國水利水電出版社,2012.

[7] KIM J K, IH J G.Prediction of sound level at high-frequency bands by means of a simplified boundary element method[J].Journal of the Acoustical Society of America, 2002,112(6): 2645-2655.

Principle to select element size of FE-BE algorithm for a ring stiffened cylindrical shell

HUANG Zhenwei, ZHOU Qidou, JI Gang, MENG Qingchang, DUAN Jiaxi

(Department of Ship and Engineering, Naval Univ.of Engineering, Wuhan 430033, China)

Investigating the principle to select element size of structural wet surface in the FE-BE algorithm has a great significance.The structural bending wavelength of the main vibration mode of a finite length ring stiffened cylinder was computed with the wave number spectrum method to provide some principle for selection of element size of a structural wet surface.Numerical computation showed the effectiveness of the principle.

ring stiffened cylindrical shell; FE-BE algorithm; vibration and sound radiation; element size

國家自然科學基金(51309230)

2015-10-27 修改稿收到日期：2016-02-19

黃振衛 男，博士生，1986年生

周其斗 男，教授，博士生導師，1962年生

TH1

A

10.13465/j.cnki.jvs.2017.05.023