郭文杰1)3) 李天勻1)2)3) 朱翔1)2)3) 屈凱旸4)

1)(華中科技大學(xué)船舶與海洋工程學(xué)院,武漢 430074)

2)(高新船舶與深海開(kāi)發(fā)裝備協(xié)同創(chuàng)新中心,上海 200240)

3)(船舶與海洋水動(dòng)力湖北省重點(diǎn)實(shí)驗(yàn)室,武漢 430074)

4)(上海交通大學(xué)船舶海洋與建筑工程學(xué)院,上海 200240)

(2017年12月18日收到;2018年2月7日收到修改稿)

1 引 言

圓柱殼結(jié)構(gòu)以其優(yōu)異的幾何特性、力學(xué)特性被廣泛應(yīng)用于土木、化工、航空航天、海洋工程等諸多領(lǐng)域,比如輸油管道、儲(chǔ)液罐、潛艇耐壓殼等.20世紀(jì),國(guó)內(nèi)外針對(duì)圓柱殼-流場(chǎng)耦合系統(tǒng)聲振問(wèn)題已有大量的文獻(xiàn)報(bào)道[1?4],也形成了較為完善的理論體系,在這些工作中,大多數(shù)的研究均針對(duì)比如全充液或者全浸沒(méi)問(wèn)題(無(wú)限流域)這種殼體與流場(chǎng)完全耦合的振動(dòng)、聲學(xué)問(wèn)題的研究.但是對(duì)于實(shí)際的工程問(wèn)題,流域往往是存在邊界的,尤其是當(dāng)結(jié)構(gòu)靠近流體邊界時(shí),或者結(jié)構(gòu)內(nèi)帶有自由液面時(shí),如果忽略流體邊界的影響,計(jì)算結(jié)果往往是有偏差的[5].

自由液面作為一類(lèi)典型的邊界條件,其與圓柱殼結(jié)構(gòu)的耦合動(dòng)力學(xué)特性吸引了許多學(xué)者對(duì)其開(kāi)展研究.Huang[6]基于鏡像原理研究了平面波入射下二維圓柱殼的散射聲場(chǎng),Hasheminejad和Azarpeyvand[7]進(jìn)一步將該方法延展到研究平面波入射下1/4空間二維圓柱的聲散射問(wèn)題.白振國(guó)等[8]采用鏡像原理建立了有限水深環(huán)境中二維圓柱殼的振動(dòng)聲學(xué)物理模型,并探討了潛深對(duì)聲場(chǎng)分布和衰減特性的影響規(guī)律.郭文杰等[9,10]基于鏡像原理和Graf加法定理建立了有限水深下有限長(zhǎng)圓柱殼-流體耦合振動(dòng)的解析模型,并基于傅里葉變換及穩(wěn)相法提出了該系統(tǒng)遠(yuǎn)場(chǎng)聲壓的快速準(zhǔn)確預(yù)報(bào)方法.

上述研究主要針對(duì)殼體結(jié)構(gòu)完全浸沒(méi)在水中的情況,這一類(lèi)問(wèn)題中,結(jié)構(gòu)完全與流場(chǎng)耦合,其聲振耦合方程的描述相對(duì)容易.實(shí)際上,還有一類(lèi)問(wèn)題,就是結(jié)構(gòu)與流場(chǎng)部分耦合的情況,比如系泊狀態(tài)下的潛艇、艦船等.對(duì)于系泊狀態(tài)下圓柱殼聲固耦合分析的研究工作則相對(duì)較少,尤其是圓柱殼軸線平行于自由液面時(shí)的解析或半解析研究更是鮮有報(bào)道,這是因?yàn)橄挡礌顟B(tài)下流體載荷解析表達(dá)式和耦合方程更難以得到.為了解決殼體部分充液的振動(dòng)問(wèn)題,Amabili[11]給出了兩種近似的方法,一種是用以圓柱圓心為原點(diǎn)構(gòu)成的扇形邊界近似替代自由液面的方法,但是這種方法僅適用于浸沒(méi)角度較小的情況;另一種是用部分環(huán)狀區(qū)域代替原始邊界,但這種方法也僅適用于浸沒(méi)角度小于π的工況.值得注意的是,Amabili提出的第一類(lèi)方法是可以推廣到殼體部分浸沒(méi)問(wèn)題中[12],并明確指出浸沒(méi)角度的適用范圍僅為?π/8—π/8;但是第二類(lèi)方法無(wú)法推廣到外流場(chǎng).基于Amabili提出的自由液面近似處理方法,Ye等[13]研究了系泊狀態(tài)下無(wú)限長(zhǎng)圓柱殼聲振特性.

此外,其他一些學(xué)者也研究過(guò)部分浸沒(méi)或者部分充液的殼體耦合振動(dòng)問(wèn)題,Selmane和Lakis[14]將圓柱殼沿周向分割成微段,再結(jié)合波動(dòng)法研究了部分充液時(shí)其流固耦合特性.Ergin和Temarel[15]基于瑞利-李茲法研究了部分充液或者部分浸沒(méi)時(shí)圓柱殼固有頻率.Krishna和Ganesan[16]基于多項(xiàng)式級(jí)數(shù)及有限元法計(jì)算了部分充液圓柱殼固有頻率以及附連水質(zhì)量.但是這些研究工作實(shí)際上僅僅考慮了濕表面流體的影響而并未考慮自由液面邊界條件的聲學(xué)效應(yīng),尚無(wú)法揭示自由液面對(duì)聲振特性的影響.另外還有一些學(xué)者通過(guò)實(shí)驗(yàn)的手段分析了部分充液的圓管自振特性[17]或者附連水質(zhì)量效應(yīng)[18],但依然無(wú)法從機(jī)理上求解這類(lèi)問(wèn)題.

需要強(qiáng)調(diào)的是,在殼體與流場(chǎng)部分耦合的問(wèn)題中,半充液或者半浸沒(méi)情況是個(gè)特例,因?yàn)檫@種模型下自由液面剛好在殼體橫截面的水平坐標(biāo)軸上,容易得到流體載荷的解析表達(dá)式,可采用正弦三角級(jí)數(shù)[19]來(lái)自動(dòng)滿足自由液面的邊界條件.基于這類(lèi)解析表達(dá)式,Li等[20]對(duì)半浸狀態(tài)下全充液圓柱殼的聲輻射特性進(jìn)行分析,李天勻等[21]針對(duì)有限長(zhǎng)半充液圓柱殼的自由振動(dòng)及受迫振動(dòng)特性開(kāi)展研究.盡管本文工作是要解決部分浸沒(méi)這類(lèi)更為普遍的工況,而文獻(xiàn)[19—21]中的研究工作僅能解決半浸沒(méi)這種特殊工況.但是,文獻(xiàn)[19—21]中關(guān)于在半浸沒(méi)情況下通過(guò)采用正弦三角級(jí)數(shù)來(lái)自動(dòng)滿足自由液面下的聲壓釋放條件的思路卻可以借鑒到本文的研究工作中來(lái).

針對(duì)更一般的殼體與流場(chǎng)部分耦合,如部分浸沒(méi)等問(wèn)題,此時(shí)自由液面與圓柱殼的圓心不共面.這給問(wèn)題的描述和求解帶來(lái)了很大的挑戰(zhàn).為解決這個(gè)問(wèn)題,本文將聲壓函數(shù)和殼體位移函數(shù)建立在不同的坐標(biāo)系下.具體來(lái)講,將聲場(chǎng)坐標(biāo)系的原點(diǎn)建立在自由液面上,并采用正弦三角級(jí)數(shù)來(lái)滿足自由液面上的聲壓釋放邊界條件;將殼體運(yùn)動(dòng)方程建立在以圓柱圓心為坐標(biāo)原點(diǎn)的坐標(biāo)系下.由此,便可以分別得到聲壓函數(shù)與殼體位移函數(shù)在各自坐標(biāo)系下的解析表達(dá)式,為解析求解部分浸沒(méi)問(wèn)題奠定了基礎(chǔ).然后,再利用Galerkin法處理結(jié)構(gòu)與聲場(chǎng)部分聲固耦合的界面上的速度連續(xù)條件,并通過(guò)兩坐標(biāo)間的幾何關(guān)系求得聲壓幅值與殼體位移幅值之間的關(guān)系矩陣,最終便可以求解該聲固耦合系統(tǒng)方程.本文方法的提出可以便捷有效地解決系泊狀態(tài)下圓柱殼的聲振問(wèn)題,豐富了該問(wèn)題理論研究的內(nèi)涵,也為解析求解浮態(tài)問(wèn)題提供了新的思路.

2 理論分析

2.1 物理模型

為了便于研究,本文假設(shè)圓柱殼軸向是無(wú)限長(zhǎng)的,且激勵(lì)力沿軸向是均勻分布的,由此本文的數(shù)學(xué)物理模型是一個(gè)典型的平面應(yīng)變模型(即二維模型).二維圓柱殼厚度為h,中面半徑為Rs,楊氏模量為E,泊松比為μ;密度為ρs,部分浸沒(méi)于密度為ρf;聲速為cf的流體中,殼體軸系與自由液面平行.如圖1所示,以殼體圓心O為坐標(biāo)原點(diǎn)的極坐標(biāo)系(r,φ)為結(jié)構(gòu)坐標(biāo)系,φ的取值范圍為0—2π,以殼體圓心正上方與自由液面的交點(diǎn)Q為原點(diǎn)的極坐標(biāo)系(R,θ)為聲場(chǎng)坐標(biāo)系,θ的取值范圍為0—π,

圖1 物理模型坐標(biāo)圖Fig.1.Coordinate figure of the physical model.

定義Q點(diǎn)到O點(diǎn)距離為浸沒(méi)深度H,浸沒(méi)角度α滿足sinα=H/Rs,當(dāng)自由液面在殼體圓心下方時(shí),H取值為負(fù)值.兩類(lèi)坐標(biāo)系與平面上任一點(diǎn)的夾角定義為β.當(dāng)自由液面在殼體圓心上方時(shí)β取值為正值,在下方時(shí)β取值為負(fù)值.徑向激勵(lì)力幅值為F0,激勵(lì)角度為φ0.

2.2 聲學(xué)邊界條件

本文的物理模型是典型的聲固耦合模型,而這項(xiàng)研究工作的難點(diǎn)和重點(diǎn)之一就是得到滿足所有對(duì)應(yīng)的聲學(xué)邊界條件的聲壓解析表達(dá)式.

首先,聲壓p必須滿足聲學(xué)Helmholtz方程:

其中kf=ω/cf為聲波波數(shù),ω=2πf為角頻率,f表示頻率,?2表示拉普拉斯算子.

其次,聲壓表達(dá)式還要滿足無(wú)限遠(yuǎn)處Sommerfeld輻射條件:

最后,自由液面處的聲壓需滿足聲壓釋放條件(因?yàn)檠芯康念l率相對(duì)較高,自由表面的重力波可以忽略不計(jì);此外由于自由液面以上空氣密度遠(yuǎn)小于水,故作為真空處理)

實(shí)際上,當(dāng)聲場(chǎng)坐標(biāo)系原點(diǎn)建立在自由液面上時(shí),可以通過(guò)采用正弦三角級(jí)數(shù)來(lái)自動(dòng)滿足自由液面聲壓為零的邊界條件[19],具體形式為

其中m為正弦三角級(jí)數(shù)的序數(shù),Pm(R)為對(duì)應(yīng)的聲壓幅值函數(shù).

對(duì)于自由液面上的任意點(diǎn),角度θ=0或π,將其代入(4)式,可以得到自由液面聲壓的表達(dá)式:

顯然此類(lèi)正弦三角級(jí)數(shù)可以用來(lái)滿足自由液面聲壓釋放條件,此外采用三角函數(shù)也更有利于應(yīng)用分離變量法來(lái)求解Helmholtz方程.由此可以得到聲壓幅值函數(shù)Pm(R)的解析表達(dá)式:

其中H(1)m()為第m階第一類(lèi)漢克爾函數(shù),Am為聲壓幅值.

由于第一類(lèi)漢克爾函數(shù)在遠(yuǎn)場(chǎng)自動(dòng)滿足Sommerfeld輻射條件,因此將(6)式代入(4)式可以得到滿足以上聲學(xué)邊界條件的聲壓解析表達(dá)式:

2.3 殼體運(yùn)動(dòng)方程

得到聲壓的解析表達(dá)式后,接下來(lái)需要建立殼體運(yùn)動(dòng)方程.本節(jié)采用二維Flügge薄殼理論[22],具體的方程為(為簡(jiǎn)潔,后文略去簡(jiǎn)諧時(shí)間項(xiàng)exp(?iωt)):

其中v和w分別是殼體中面切向和徑向位移,fp表示作用在圓柱殼表面的聲載荷,f0表示外激勵(lì)力載荷,[L]是二維Flügge薄殼方程中的微分算子矩陣,具體如下:

由于殼體和聲介質(zhì)是部分耦合,因此(8)式中fp應(yīng)表示為分段函數(shù)的形式:

假設(shè)作用于圓柱殼的激勵(lì)力是一個(gè)沿軸向均勻分布的無(wú)限長(zhǎng)徑向線力,激勵(lì)力位于結(jié)構(gòu)坐標(biāo)系的(Rs,φ0)處,則激勵(lì)力載荷可以表示為

其中δ()表示狄利克雷函數(shù).

對(duì)于圓柱殼結(jié)構(gòu),由于周向的周期性,其位移及載荷函數(shù)可以在周向展開(kāi)為傅里葉級(jí)數(shù)的形式[8]:

其中Vn和Wn分別是周向和徑向位移幅值,fpn和f0n分別表示殼體表面聲載荷fp以及激勵(lì)力載荷f0的幅值,n是傅里葉周向展開(kāi)序列數(shù).

因?yàn)?9)和(13)式是殼體表面聲載荷fp的不同形式的表達(dá)式,利用正交化處理,可得fpn的表達(dá)式:

由圖1可知,空間中任意點(diǎn)的結(jié)構(gòu)坐標(biāo)與聲學(xué)坐標(biāo)有如下關(guān)系:

因此對(duì)于流固耦合面上任一點(diǎn),令r=Rs,則(R,θ)在聲學(xué)坐標(biāo)系下的坐標(biāo)可以由(16)式解出:

將(17)式代入(15)式即可對(duì)fnp進(jìn)行求解,但是由于積分中包含漢克爾函數(shù),無(wú)法直接進(jìn)行積分計(jì)算,因此本文采用離散求和的形式來(lái)近似計(jì)算.首先將積分域均分為K段,然后取各段中點(diǎn)值代入求和公式中:

其中F(φ)表示需要進(jìn)行積分的函數(shù),通過(guò)大量算例表明K取100時(shí)收斂性已經(jīng)很好.

同理,(10)和(14)式是外激勵(lì)力f0的不同形式的表達(dá)式,利用正交化處理,可以得到f0n的表達(dá)式:

然后將(11)—(14)式代入到(8)式中,并進(jìn)行正交化處理,可以得到解耦后的殼體運(yùn)動(dòng)方程:

其中矩陣[T]的元素如下:T11=?2?n2,T12=T21=in,T22=1+K+Kn4?2Kn2??2.是無(wú)量綱頻率.

由(20)式可以得到僅與徑向位移幅值相關(guān)的控制方程:

其中In=T11/det(T),In與n有關(guān),det(T)表示矩陣[T]的行列式.

很明顯,求解控制方程(21)的關(guān)鍵在于得到徑向位移幅值Wn和殼體表面聲載荷幅值fpn之間的關(guān)系.因此需要通過(guò)殼體表面流體與結(jié)構(gòu)速度連續(xù)條件來(lái)求解該問(wèn)題,在聲固耦合交界面處.具體方程為[10]:

因?yàn)楹茈y直接求解速度連續(xù)方程,采用Galerkin法進(jìn)行處理,可選擇的權(quán)函數(shù)有兩類(lèi),一類(lèi)是殼體徑向位移的周向函數(shù),另一類(lèi)是聲壓的周向函數(shù):

因此方程(22)可以轉(zhuǎn)變?yōu)镚alerkin積分的弱形式:

其中N是截?cái)囗?xiàng)數(shù),即需要構(gòu)造2N+1個(gè)積分方程.

由(24)式可以得到徑向位移幅值和聲壓幅值之間的關(guān)系:

其 中[Vs]和[Vp]均 是2N+1階 方 陣,{Wn}和{Am}分別表示徑向位移幅值向量和聲壓幅值向量,且{Wn}=[W?N,W?N+1,···,WN?1,WN]T,{Am}=[A1,A2,···,A2N,A2N+1]T,上標(biāo)T表示轉(zhuǎn)置.

根據(jù)(16)式中兩類(lèi)坐標(biāo)系的幾何關(guān)系,可以將聲壓沿徑向的導(dǎo)數(shù)轉(zhuǎn)換到聲學(xué)坐標(biāo)系下:

其中夾角β=3π/2?θ?φ.

由此,當(dāng)權(quán)函數(shù)選擇殼體徑向位移的周向函數(shù)exp(inφ)時(shí),[Vs]和[Vp]中每一個(gè)元素的具體表達(dá)式如下:

其中a,b分別表示矩陣的行和列.

當(dāng)權(quán)函數(shù)選擇聲壓的周向函數(shù)sin(mθ)時(shí),[Vs]和[Vp]中每一個(gè)元素的具體表達(dá)式如下:

由于(27)或(28)式中積分包含漢克爾函數(shù),無(wú)法直接進(jìn)行積分計(jì)算,因此也采用(18)式中離散求和的方法來(lái)近似計(jì)算.

為了方便求解耦合方程,將(21)式改寫(xiě)為矩陣運(yùn)算方程:

其中{f0n}=F0/2π{exp(iNφ0),exp[i(N?1)φ0],···,exp(?iNφ0)}T;[G]矩陣為對(duì)角矩陣,[G]j,j=Ij?1?N;{fpn}={fp,?N,fp,?N+1,···,fp,N}T.

同時(shí),(15)式中{fpn}也可表示為矩陣的形式:{fpn}=[Tp]{Am}.(30)其中矩陣[Tp]中每一個(gè)元素的具體表達(dá)式如下:

將(25)和(30)式代入方程(29)中,可以求解出徑向位移幅值{Wn}:

其中矩陣[J]是2N+1階單位矩陣,并且在得到徑向位移幅值{Wn}后,根據(jù)(25)式可以求解出{Am},從而可以計(jì)算聲壓.

前文求解受迫振動(dòng)時(shí),是已知激勵(lì)力及激勵(lì)頻率,求響應(yīng)(徑向位移幅值).當(dāng)求解自由振動(dòng)時(shí),并沒(méi)有激勵(lì)源,固有頻率是要求解的未知量.由此(32)式可以表示為

其中{0}表示零向量,顯然(33)式是個(gè)典型的特征值問(wèn)題,即求解角頻率ω.

具體來(lái)講就是通過(guò)定義(33)式中系數(shù)矩陣行列式值為零來(lái)求解角頻率:

其中F(ω)表示系數(shù)矩陣的行列式值.

因?yàn)榫仃嘯J],[Vp],[Tp]中均含有角頻率ω,實(shí)際上矩陣的行列式值F(ω)=0是個(gè)關(guān)于角頻率隱式表達(dá)的超越方程,難以直接解出ω.因此本文采用搜根的方法進(jìn)行求解,通過(guò)設(shè)置合適的搜根步長(zhǎng),逐步增大輸入的角頻率,當(dāng)行列式值過(guò)零點(diǎn)時(shí),輸出對(duì)應(yīng)的角頻率,這個(gè)角頻率即為系統(tǒng)的固有角頻率.由此可以將超越方程的各階角頻率從小到大逐一求出,再根據(jù)可以得到各階固有頻率.

3 數(shù)值計(jì)算

模型參數(shù):半徑Rs=0.18 m,厚度h=0.001 m,殼體密度ρ=7850 kg/m3,楊氏模量E=206 GPa,泊松比μ=0.3,流體密度ρf=1025 kg/m3,流體聲速cf=1500 m/s.

3.1 收斂性分析

為了說(shuō)明方法的收斂性,取無(wú)量綱浸沒(méi)深度H/Rs=?0.5,H/Rs=0和H/Rs=0.5,分別計(jì)算在激勵(lì)力頻率200,400,800 Hz時(shí)徑向均方根振速Vm隨截?cái)鄶?shù)N的變化規(guī)律,定義

其中vn=?w/?t表示徑向速度;則均方根速度級(jí) VML=20lg(Vm/V0),其中基準(zhǔn)速度V0=10?6m/s. 激勵(lì)力幅值F0=1 N,激勵(lì)角度φ0= π/4. 定義復(fù)楊氏模量E′=E(1+iη),結(jié)構(gòu)阻尼η=0.01.計(jì)算中權(quán)函數(shù)選擇位移的周向展開(kāi)函數(shù).

圖2 均方根速度級(jí)VML收斂性分析Fig.2.Convergence analysis of the root mean square velocity levels.

從圖2可見(jiàn),均方根速度級(jí)隨著截?cái)囗?xiàng)數(shù)N的增大很快趨于穩(wěn)定;并且頻率越高,達(dá)到收斂時(shí)截?cái)囗?xiàng)數(shù)N的取值越大.由圖2可知,對(duì)于頻率小于800 Hz時(shí)的分析計(jì)算,N取16時(shí)已足夠收斂.

3.2 本文方法的適用性分析

對(duì)于Amabili[12]提出的自由液面近似處理方法,浸沒(méi)角度α的適用范圍僅為?π/8到π/8(無(wú)量綱浸沒(méi)深度H/Rs約為?0.38—0.38,無(wú)量綱深度為負(fù)值表示殼體的圓心在自由液面上方),而本文方法的主要優(yōu)勢(shì)就是模型中殼體部分浸沒(méi)時(shí)的浸沒(méi)深度范圍更廣.為了驗(yàn)證本文的方法,無(wú)量綱浸沒(méi)深度取值為?0.9—0.9,分別計(jì)算兩種不同權(quán)函數(shù)下首階固有頻率值,并與有限元軟件Comsol仿真計(jì)算結(jié)果進(jìn)行對(duì)比.其中有限元模型如圖3所示,流域以聲學(xué)坐標(biāo)原點(diǎn)為中心,半徑取2 m,用完美匹配層模擬無(wú)限遠(yuǎn)聲學(xué)邊界,匹配層厚度取0.05 m.網(wǎng)格包含7980個(gè)域單元和756個(gè)邊界單元.

圖3 有限元模型Fig.3.The finite element model.

表1 不同浸沒(méi)深度下首階固有頻率對(duì)比(單位為Hz)Table 1.Comparison of the fundamental frequency with dif f erent immerged depth(in Hz).

從表1可以看出,選擇位移函數(shù)或者選擇聲壓函數(shù)作為權(quán)函數(shù)時(shí),首階固有頻率計(jì)算結(jié)果符合很好.并且無(wú)論選擇哪類(lèi)權(quán)函數(shù),無(wú)量綱浸沒(méi)深度H/Rs均可以從?0.9變化到0.8(當(dāng)無(wú)量綱浸沒(méi)深度過(guò)大或者過(guò)小時(shí)計(jì)算結(jié)果難以收斂),這表明本文的浸深適用范圍非常廣,方法更具有一般性.

此外,從表1還可以看出,本文方法計(jì)算結(jié)果與Comsol仿真計(jì)算結(jié)果符合良好,而且不同方法下首階固有頻率隨浸沒(méi)深度的變化規(guī)律也是一致的.這是由于隨著殼體浸沒(méi)深度的增大,流固耦合面增大,附連水質(zhì)量也相應(yīng)增加,因此固有頻率會(huì)逐漸減小.

3.3 準(zhǔn)確性驗(yàn)證

3.3.1 自由振動(dòng)的驗(yàn)證

為說(shuō)明本文方法計(jì)算自由振動(dòng)問(wèn)題的準(zhǔn)確性,分別取無(wú)量綱浸沒(méi)深度H/Rs=?0.7和H/Rs=0.7,計(jì)算系統(tǒng)前十階固有頻率,并與有限元軟件Comsol仿真計(jì)算結(jié)果進(jìn)行對(duì)比.定義固有頻率的相對(duì)誤差Error=|f1?f2|/f2×100%.

表2 H/Rs=?0.7時(shí)前十階固有頻率對(duì)比(單位為Hz)Table 2.Comparison of the the first ten order natural frequencies when H/Rs=?0.7(in Hz).

表3 H/Rs=0.7時(shí)前十階固有頻率對(duì)比(單位為Hz)Table 3.Comparison of the the first ten order natural frequencies when H/Rs=0.7(in Hz).

從表2和表3可以看出,本文方法計(jì)算得到的前10階固有頻率值與Comsol仿真計(jì)算結(jié)果符合良好,最大相對(duì)誤差不超過(guò)1%,說(shuō)明采用本文方法計(jì)算固有頻率是準(zhǔn)確可靠的.

3.3.2 受迫振動(dòng)的準(zhǔn)確性驗(yàn)證

在分析完自由振動(dòng)的準(zhǔn)確性之后,進(jìn)一步分析受迫振動(dòng)的準(zhǔn)確性. 取無(wú)量綱浸沒(méi)深度H/Rs=?0.5和H/Rs=0.5,計(jì)算在激勵(lì)力頻率50—500 Hz時(shí)測(cè)點(diǎn)徑向速度(頻率間隔10 Hz),定義徑向速度級(jí)VL=20×log(|V|/V0),其中V為徑向速度,基準(zhǔn)速度V0=10?6m/s.激勵(lì)力幅值F0=1 N,激勵(lì)角度φ0=0.測(cè)點(diǎn)位于周向角φ=π處.

從圖4可以看出,本文計(jì)算結(jié)果和Comsol仿真計(jì)算結(jié)果整體符合良好,說(shuō)明本文方法分析受迫振動(dòng)是準(zhǔn)確可靠的.但是當(dāng)頻率較高時(shí)(以圖4為例,大于400 Hz),誤差逐漸增大.主要的原因可能是有限元計(jì)算聲固耦合問(wèn)題時(shí),隨著頻率增大,對(duì)網(wǎng)格密度的要求也提高,計(jì)算精度會(huì)降低.

圖4 不同方法下徑向速度級(jí)對(duì)比 (a)H/Rs=?0.5;(b)H/Rs=0.5Fig.4.Comparison of the radial velocity levels with dif f erent methods:(a)H/Rs=?0.5;(b)H/Rs=0.5.

3.3.3 聲場(chǎng)求解的準(zhǔn)確性驗(yàn)證

為進(jìn)一步說(shuō)明本文方法計(jì)算聲場(chǎng)也是準(zhǔn)確的,取無(wú)量綱浸沒(méi)深度H/Rs=?0.4和H/Rs=0.4,計(jì)算激勵(lì)頻率為50 Hz時(shí)聲壓幅值云圖,并與Comsol仿真計(jì)算結(jié)果進(jìn)行對(duì)比.其中激勵(lì)力幅值F0=1 N,激勵(lì)角度φ0=0,云圖尺寸為1.2 m×0.6 m.

從圖5和圖6可以看出,本文方法計(jì)算得到的聲壓云圖和有限元軟件Comsol仿真計(jì)算結(jié)果符合良好,由此可以說(shuō)明本文方法計(jì)算聲場(chǎng)是準(zhǔn)確可靠的.

圖5 H/Rs=?0.4時(shí)不同方法下聲壓幅值云圖對(duì)比 (a)本文方法;(b)Comsol仿真Fig.5.Comparison of the sound pressure contour map with dif f erent methods when H/Rs=?0.4:(a)Present method;(b)Comsol.

圖6 H/Rs=0.4時(shí)不同方法下聲壓幅值云圖對(duì)比 (a)本文方法;(b)Comsol仿真Fig.6.Comparison of the sound pressure contour map with dif f erent methods when H/Rs=0.4:(a)Present method;(b)Comsol.

此外,值得一提的是,本文方法計(jì)算效率也非常高,以圖5或圖6中聲壓云圖的計(jì)算為例,在Matlab中僅僅需要不到2 s即可計(jì)算出精確穩(wěn)定的結(jié)果.

4 討 論

在驗(yàn)證了本文方法計(jì)算聲固耦合系統(tǒng)振動(dòng)及聲壓求解均準(zhǔn)確可靠之后,進(jìn)一步分析與討論部分浸沒(méi)殼體的模態(tài)振型、無(wú)量綱浸沒(méi)深度對(duì)受迫振動(dòng)的影響以及遠(yuǎn)場(chǎng)聲壓的指向性特征.

4.1 模態(tài)振型分析

為了更直觀地揭示聲固部分耦合系統(tǒng)自振特性,取無(wú)量綱浸沒(méi)深度H/Rs=0.7和無(wú)限域時(shí)計(jì)算得到的前4階模態(tài)振型及固有頻率.

對(duì)比表4和表5可以看出,部分浸沒(méi)工況下振型與無(wú)限域振型有差異.以H/Rs=0.7時(shí)第一階振型為例,振型函數(shù)約為cos2φ+0.25·cosφ?0.05·cos3φ,而對(duì)于無(wú)限域殼體,對(duì)應(yīng)的振型函數(shù)為cos2φ(規(guī)則的周向波型).這是由于流體在周向分布不均,殼體與結(jié)構(gòu)的部分耦合破壞了圓柱殼的周向的對(duì)稱性,故規(guī)則的周向波之間會(huì)發(fā)生互耦,形成復(fù)雜的振型函數(shù).

表4 H/Rs=0.7時(shí)前4階模態(tài)振型Table 4.Modal shapes of the first four orders when H/Rs=0.7.

表5 無(wú)限域時(shí)前4階模態(tài)振型Table 5.Modal shapes of the first four orders when in infinite fluid.

此外,由于自由液面的存在,系統(tǒng)僅有惟一對(duì)稱軸,對(duì)稱和反對(duì)稱模態(tài)固有頻率也存在差異.以H/Rs=0.7算例下前兩階固有頻率為例,分別為7.91 Hz和8.86 Hz,有明顯差異;但是對(duì)于無(wú)限域情況,由于系統(tǒng)具有周向?qū)ΨQ性,對(duì)稱和反對(duì)稱模態(tài)固有頻率是相同的,二者存在明顯區(qū)別.

另外,部分浸沒(méi)工況下,對(duì)稱和反對(duì)稱模態(tài)振型函數(shù)之間的周向波互耦程度也不盡相同.仍以H/Rs=0.7工況為例,第二階振型函數(shù)(反對(duì)稱)約為sin2φ+0.38·sinφ?0.08·sin3φ,函數(shù)各成分之間的比值由對(duì)稱模態(tài)時(shí)1:0.25:?0.05變?yōu)榉磳?duì)稱時(shí)1:0.38:?0.08.而對(duì)于無(wú)限域情況,對(duì)應(yīng)的振型函數(shù)為sin2φ,這也是系泊狀態(tài)下自振特性區(qū)別于無(wú)限域工況的又一特征.

4.2 均方根振速頻譜分析

3.2節(jié)中提到,隨著浸沒(méi)深度增大,流固耦合面會(huì)增大,附連水質(zhì)量也會(huì)相應(yīng)增加.為進(jìn)一步解釋這個(gè)觀點(diǎn),本文選取了H/Rs=?0.6,?0.3,0,0.3,0.6這五個(gè)深度,對(duì)比分析各工況下均方根速度級(jí)的頻譜曲線.激勵(lì)力幅值F0=1 N,激勵(lì)角度φ0=0,激勵(lì)頻率為1—100 Hz,掃頻間隔1 Hz.

從圖7可以看出,隨著浸沒(méi)深度增大(液面升高),共振峰均向低頻移動(dòng).這也是因?yàn)榻](méi)深度增大使得附連水質(zhì)量增大,從而增加了系統(tǒng)的總質(zhì)量,導(dǎo)致共振頻率降低,頻譜曲線整體左移.

4.3 遠(yuǎn)場(chǎng)聲壓指向性分析

系泊狀態(tài)下的艦船等目標(biāo)的水下輻射噪聲對(duì)于其隱蔽性有著重要意義.因此,本節(jié)進(jìn)一步開(kāi)展二維系泊圓柱殼遠(yuǎn)場(chǎng)聲壓指向性特征的研究.

取無(wú)量綱浸沒(méi)深度H/Rs=?0.5,0,0.5,分別計(jì)算激勵(lì)頻率f=100,200 Hz時(shí)遠(yuǎn)場(chǎng)聲壓級(jí)指向性圖.其中激勵(lì)力幅值F0=1 N,激勵(lì)位置φ0=π/4.遠(yuǎn)場(chǎng)點(diǎn)取自聲學(xué)坐標(biāo)系,半徑R=1000 m,角度θ取0—π,取值間隔θ/180.定義聲壓級(jí)SPL=20lg(|p|/p0),其中基準(zhǔn)聲壓p0=10?6m/s.

從圖8可以看出,所有聲壓指向性曲線都在角度θ=π/2時(shí)(正下方)取最大值.而且盡管激勵(lì)位置并不在殼體的垂直對(duì)稱軸上,但聲壓級(jí)關(guān)于垂直對(duì)稱軸呈對(duì)稱分布.以下對(duì)這些現(xiàn)象進(jìn)行物理解釋.

圖7 不同浸沒(méi)深度下均方根速度級(jí)頻譜曲線 (a)H/Rs=?0.6,?0.3,0;(b)H/Rs=0,0.3,0.6Fig.7.Spectrum curves of the root mean square velocity levels at dif f erent immersion depths:(a)H/Rs=?0.6,?0.3,0;(b)H/Rs=0,0.3,0.6.

圖8 不同浸沒(méi)深度下聲壓級(jí)指向性 (a)f=100 Hz;(b)f=200 HzFig.8.Directivity of the sound pressure levels at dif f erent immersion depths:(a)f=100 Hz;(b)f=200 Hz.

因?yàn)檠芯康膱?chǎng)點(diǎn)位于遠(yuǎn)場(chǎng),從幾何上講,殼體濕表面上任意一點(diǎn)到場(chǎng)點(diǎn)的距離可以認(rèn)為是近似相同的.由此水下的輻射面可以等效為一個(gè)點(diǎn)源,而且這個(gè)點(diǎn)源距離自由液面的距離在0—2Rs之間(這個(gè)距離與速度分布有關(guān),是未知的).

基于鏡像原理[8],自由液面對(duì)水下聲場(chǎng)的作用可以通過(guò)在自由液面的另一側(cè)構(gòu)造等距的(虛源到液面距離)以及反相位的虛源來(lái)實(shí)現(xiàn).則水下聲場(chǎng)可以認(rèn)為是由實(shí)源和虛源共同貢獻(xiàn),且虛源和實(shí)源共同作用下聲壓在自由液面上滿足為零的條件.實(shí)際上,由于實(shí)源和虛源距離較近,在低頻下,這樣的物理模型構(gòu)成了偶極子模型.由此,可以得到聲偶極模型的數(shù)學(xué)表達(dá)式[8]:

其中A表示點(diǎn)源的聲壓幅值,R代表場(chǎng)點(diǎn)到聲學(xué)坐標(biāo)原點(diǎn)的距離,D表示點(diǎn)源與鏡像源距離(D/Rs<4).

雖然D值是未知的,但并不影響分析聲壓的指向性.因?yàn)楫?dāng)激勵(lì)頻率不太高時(shí),kf和D的乘積是小于π/2的.因此當(dāng)遠(yuǎn)場(chǎng)點(diǎn)位于殼體截面圓心正下方時(shí)(θπ/2),從(35)式可以看出,聲壓幅值將取最大值.此外,從(35)式中正弦三角函數(shù)的數(shù)學(xué)性質(zhì)也可以看出,聲壓級(jí)曲線關(guān)于對(duì)稱軸(θ=π/2)將呈對(duì)稱分布,這也很好地契合了圖8中聲壓級(jí)指向性曲線的特征.

5 結(jié) 論

提出了一種求解系泊狀態(tài)下圓柱殼聲振問(wèn)題的新方法.通過(guò)將聲壓與殼體位移建立在不同的坐標(biāo)系下,方便地得到了其解析表達(dá)式,然后再利用Galerkin法以及坐標(biāo)變換關(guān)系處理聲固耦合界面的速度連續(xù)條件,最終可以很便捷、高效地求解該浮態(tài)系統(tǒng)的聲振特性,也為求解半空間中彈性結(jié)構(gòu)與聲場(chǎng)部分耦合的聲振問(wèn)題提供了新的思路.具體結(jié)論如下.

1)通過(guò)對(duì)部分浸沒(méi)殼體的自由、受迫振動(dòng)和輻射聲壓的求解并與有限元數(shù)值解進(jìn)行對(duì)比分析,驗(yàn)證了本文方法的準(zhǔn)確性,且本文方法能適用于浸沒(méi)深度更大范圍變化的一般情況.

2)隨著自由液面的升高和下潛深度的增大,流體與結(jié)構(gòu)耦合交界面積增大,導(dǎo)致附連水質(zhì)量相應(yīng)增大,從而整個(gè)系統(tǒng)的附連水質(zhì)量也會(huì)增大,因此同階次固有頻率會(huì)逐漸減小.

3)部分浸沒(méi)工況下自振特性分析時(shí),流體與結(jié)構(gòu)處于部分耦合的狀態(tài),從而導(dǎo)致由結(jié)構(gòu)傳遞到聲場(chǎng)中的聲波碰到聲學(xué)邊界產(chǎn)生回波,疊加形成復(fù)雜的聲場(chǎng).流體介質(zhì)對(duì)聲場(chǎng)的作用以聲負(fù)載的形式作用于結(jié)構(gòu)表面,導(dǎo)致結(jié)構(gòu)的周向波發(fā)生互耦,由此結(jié)構(gòu)的周向模態(tài)振型不再是規(guī)則的波型.

4)由于自由液面的存在使得整個(gè)耦合系統(tǒng)僅有惟一對(duì)稱軸,從而破壞了圓柱殼振動(dòng)在周向的任意對(duì)稱性.因此結(jié)構(gòu)的對(duì)稱和反對(duì)稱模態(tài)的固有頻率會(huì)有明顯差異.

5)由于類(lèi)偶極子效應(yīng),輻射聲波以及由邊界反射的聲波發(fā)生干涉.當(dāng)激勵(lì)頻率相對(duì)較低時(shí),根據(jù)相干波疊加函數(shù)的性質(zhì),遠(yuǎn)場(chǎng)聲壓級(jí)的最大值總是出現(xiàn)在結(jié)構(gòu)正下方.

[1]Williams W,Parke N G,Moran D A,Sherman C H 1964J.Acount.Soc.Am.36 2316

[2]Harari A,Sandman B E 1976J.Acoust.Soc.Am.60 117

[3]Fuller C R 1988J.Sound Vib.122 479

[4]Pan A,Fan J,Wang B,Chen Z G,Zheng G Y 2014Acta Phys.Sin.63 214301(in Chinese)[潘安,范軍,王斌,陳志剛,鄭國(guó)垠2014物理學(xué)報(bào)63 214301]

[5]Liu P,Liu S W,Li S 2013Ship Sci.Technol.36 36(in Chinese)[劉佩,劉書(shū)文,黎勝 2013艦船科學(xué)技術(shù) 36 36]

[6]Huang H 1981Wave Motion3 269

[7]Hasheminejad S M,Azarpeyvand M 2005J.Appl.Math.Mech.85 66

[8]Bai Z G,Wu W W,Zuo C K,Zhang F,Xiong C X 2014J.Ship Mech.18 178(in Chinese)[白振國(guó),吳文偉,左成魁,張峰,熊晨熙2014船舶力學(xué)18 178]

[9]Guo W J,Li T Y,Zhu X,Miao Y Y,Zhang G J 2017J.Sound Vib.393 338

[10]Guo W J,Li T Y,Zhu X,Zhang S 2017Chin.J.Ship Res.12 62(in Chinese)[郭文杰,李天勻,朱翔,張帥2017中國(guó)艦船研究12 62]

[11]Amabili M 1996J.Sound Vib.191 757

[12]Amabili M 1997J.Vib.Acoust.119 476

[13]Ye W B,Li T Y,Zhu X 2012Appl.Mech.Mater.170 2303

[14]Selmane A,Lakis A A 1997J.Sound Vib.5 111

[15]Ergin A,Temarel P 2002J.Sound Vib.254 951

[16]Krishna B V,Ganesan N 2006J.Sound Vib.291 1221

[17]Hidalgo J A,Gama A L,Moreire R M 2017J.Sound Vib.408 31

[18]Escaler X,Torre O D,Goggins J 2017J.Fluid Struct.69 252

[19]Li T Y,Wang P,Zhu X,Yang J,Ye W B 2017J.Vib.Acoust.139 041002

[20]Li H L,Wu C J,Huang X Q 2003Appl.Acoust.64 495

[21]Li T Y,Wang L,Guo W J,Yang G D,Zhu Xiang 2016Chin.J.Ship Res.11 106(in Chinese)[李天勻,王露,郭文杰,楊國(guó)棟,朱翔2016中國(guó)艦船研究11 106]

[22]Zhu X,Ye W B,Li T Y,Chen C 2013Ocean Eng.58 22