劉 偉, 趙千里
(常州機(jī)電職業(yè)技術(shù)學(xué)院 機(jī)械工程學(xué)院,江蘇 常州 213164)
隨著人類在探索土地、深海、太空等領(lǐng)域的不斷深入,對各類機(jī)械裝備的性能也提出了越來越苛刻的要求,輸流管路作為大部分整機(jī)系統(tǒng)不可或缺的組成部分,因其在制造材料以及截面形狀等方面具有多樣性,還因其便于生產(chǎn)制造、可通過組合構(gòu)造任意空間構(gòu)型等特點(diǎn),輸流管路被廣泛應(yīng)用于各種場合,如:空間站機(jī)械臂的液壓系統(tǒng),長距離石油、天然氣的傳輸系統(tǒng),城市的供熱管網(wǎng)系統(tǒng),中央空調(diào)的送風(fēng)系統(tǒng),農(nóng)田灌溉用水泵的抽水、送水系統(tǒng),化工行業(yè)液體原料的輸送系統(tǒng)等。
正因如此廣泛的應(yīng)用范圍,導(dǎo)致在正常服役期間,輸流管路難免會(huì)因內(nèi)部流體壓力不均、流速過快、存在雜質(zhì)、材料本身存在缺陷、外部存在激振源等因素而誘發(fā)管路系統(tǒng)發(fā)生機(jī)械振動(dòng)。在近半個(gè)世紀(jì)以來,借助飛速發(fā)展的自然科學(xué)技術(shù),學(xué)者們從多個(gè)方面對管路的流固耦合振動(dòng)問題開展了研究逐漸形成了兩大研究分支[1]:其一,建立問題的數(shù)學(xué)模型;其二,開發(fā)高效準(zhǔn)確的數(shù)值方法。
目前,人們在建立研究輸流管路流固耦合振動(dòng)問題的數(shù)學(xué)模型時(shí)多以梁模型為理論基礎(chǔ),如:以Paidoussis教授為首的科研團(tuán)隊(duì)[2-7]基于Euler-Bernoulli梁模型,從支承、材料屬性、空間構(gòu)型、激勵(lì)等方面出發(fā)建立了管路的振動(dòng)微分方程,并提出了許多重要的理論,如:Paidoussis耦合模態(tài)顫振、關(guān)于彎管軸線是否可伸長的三大理論等。
在數(shù)值方法方面,在近半個(gè)世紀(jì)以來,隨著計(jì)算機(jī)技術(shù)的不斷進(jìn)步,大量的數(shù)值方法也隨之發(fā)展起來,在求解管路振動(dòng)特性方面典型的方法主要包括:傳遞矩陣法[8],微分變換法[9],微分求積法[10]或其廣義形式[11],格林函數(shù)法[12]等。
前面提到的所有研究結(jié)果的共同點(diǎn)在于:管路的放置方向?yàn)樗交蜇Q直??紤]到在實(shí)際的工程應(yīng)用中,傾斜安置的管路并不少見,如:草坪的噴灑灌溉系統(tǒng)、噴泉的供水系統(tǒng)等。此外,很多場合下管路上需要添加附加質(zhì)量以控制系統(tǒng)的動(dòng)力學(xué)表現(xiàn),如:消防車的水炮、發(fā)電站冷卻塔內(nèi)的噴淋管路、含噴嘴的噴泉供水管路及含附加質(zhì)量的城市灑水車的出水管路等,然而,同時(shí)研究這兩個(gè)因素對管路振動(dòng)特性的影響的報(bào)道卻鮮有提及。
由于多數(shù)場合下管路兩端的支承沿管路軸線方向具有一定的長度,可以近似視為固定-固定的支承形式,因此,本文重點(diǎn)研究在該支承形式下含集中質(zhì)量的斜置輸流管路的流體誘發(fā)振動(dòng)特性。首先,以數(shù)學(xué)表達(dá)式表示了由傾斜角度和集中質(zhì)量引起的管路受力的變化,而后將它們引入Euler-Bernoulli型輸流直管基本的流體誘發(fā)振動(dòng)微分方程,得到了系統(tǒng)的動(dòng)力學(xué)偏微分方程;其次,基于作者提出的基于新型形函數(shù)的Galerkin法[13]對上述偏微分方程進(jìn)行離散,經(jīng)過推導(dǎo),獲得了計(jì)算系統(tǒng)固有頻率的特征方程;最后,研究了傾角和集中質(zhì)量對管路固有頻率及臨界流速的影響。
含集中質(zhì)量的斜置輸流管路的力學(xué)模型如圖1所示,其中管路總長以L表示,單位為m。
圖1 含集中質(zhì)量的斜置輸流管路的力學(xué)模型
一般而言,如果不計(jì)軸向張力、內(nèi)部壓力、管路材料的能量耗散、管路和流體間的阻尼等因素,輸流直管的橫向振動(dòng)微分方程可表示為
式中:E為彈性模量,N/m2;I為橫截面慣性矩,m4;w為截面的橫向位移,m;x為位置坐標(biāo),m;t為時(shí)間,s;mf為單位長度流體的質(zhì)量,kg/m;mp為單位長度管路的質(zhì)量,kg/m;U為橫截面內(nèi)流體的平均流速,m/s。
在式(1)中,等號左端從左至右依次表示管路的彎曲回復(fù)力,流體的離心力,流體的科氏力以及管路和流體的慣性力。該式適用于水平放置且在力學(xué)模型上與Euler-Bernoulli梁相近的管路系統(tǒng)。而當(dāng)管路傾斜放置時(shí),需要計(jì)入管路連同流體微元以及附加質(zhì)量的重力對橫向振動(dòng)的影響,為便于描述,將重力沿管路軸線方向和與之垂直的法向分解為兩個(gè)分力,則上述影響主要體現(xiàn)在:
(1)重力在軸線方向的分力引起的附加力[14],并且在法線方向產(chǎn)生分力;
(2)附加質(zhì)量引起慣性力發(fā)生變化和產(chǎn)生轉(zhuǎn)動(dòng)慣性力[15]。
因此,以上兩個(gè)因素引起的附加項(xiàng)可以表示為
(2)
式中:Γ=mf+mp+M0δ(x-x0);M0為附加質(zhì)量,kg;J0為轉(zhuǎn)動(dòng)慣量,kg·m2,其定義為集中質(zhì)量與它在轉(zhuǎn)動(dòng)過程中回轉(zhuǎn)半徑的平方的乘積;g為重力加速度,m/s2;x0為附加質(zhì)量的安置坐標(biāo),m;α為管路的傾斜角度,rad。
將式(2)代入式(1),于是圖1所示管路的控制方程可以表示為
從微分方程的角度來看,式(3)所描述的非齊次偏微分方程與對應(yīng)的齊次方程的通解一致。而且,從振動(dòng)力學(xué)的角度來看,考慮到管路中心線的法向分力僅能影響振動(dòng)的平衡位置而與系統(tǒng)的固有特性無關(guān),因此,為了求式(3)的通解,可將其改寫為
(4)
式(4)的無量綱形式為
式中:κ=λ+γδ(ξ-ξ0);ε=γδ(ξ-ξ0)/λ。且各參量的定義為
(6)
對于如圖1所示的兩端固定式輸流直管,其無量綱的邊界條件可以表示為
η(0,τ)=η′(0,τ)=η(1,τ)=η′(1,τ)=0
(7)
式(5)的解可由Galerkin法表示為
(8)
式中:N為形函數(shù)的個(gè)數(shù);φn為第n個(gè)形函數(shù);qn為與之對應(yīng)的時(shí)間相關(guān)項(xiàng)。
將式(8)代入式(5),可得
(9)
分別將各形函數(shù)與式(9)相乘并將結(jié)果關(guān)于ξ在區(qū)間[0,1]內(nèi)積分,經(jīng)過整理可得
(10)
其中,
經(jīng)過整理,式(10)可以表示為
(11)
式(11)的解可以表示為
q=q0exp(iωτ)
(12)
將式(12)代入式(11),經(jīng)過整理可得
[K+iωG-ω2M]q0=0
(13)
因此,為了得到q0的非平凡解,需滿足
|K+iωG-ω2M|=0
(14)
式(14)即為求解系統(tǒng)特征值的特征方程。由上述推到過程可知,為了最終求解式(14),需要知道形函數(shù)向量的具體形式,本文采用趙千里等研究中的形函數(shù)表達(dá)式,即
φn=anξn+1(1-ξ)2, (n=1, 2, …,N)
(15)
其中,
形函數(shù)向量則表示為
(N+2)(N+3)(N+4)ξN+1}T
(16)
將式(15)引入式(9)~式(14)便可得到系統(tǒng)的特征值。根據(jù)上述過程,可知特征值ωj(j=1,2,3,…)僅與E,I,L,U,mf,mp,M,x0和α等系統(tǒng)參數(shù)有關(guān)且必然為復(fù)數(shù),據(jù)孫志禮等所述,它的實(shí)部表示系統(tǒng)的固有頻率,虛部與系統(tǒng)的阻尼有關(guān)。
以某段兩端固定式輸流直管為例,其中,管路部分的已知參數(shù)包括:長度L=1 000 mm,彈性模量E=2.06×1011Pa,密度橫截面內(nèi)、外圓半徑分別為r1=12 mm,r2=14 mm,密度ρp=7 800 kg/m3,而管路內(nèi)部的流體為水,密度ρf=1 000 kg/m3,經(jīng)式(6)參量的定義計(jì)算,質(zhì)量比β=0.262。趙千里等驗(yàn)證了當(dāng)函數(shù)個(gè)數(shù)N=8時(shí),本文方法與格林函數(shù)法的解已經(jīng)十分相近,隨著N的增加,本方法的計(jì)算結(jié)果將無限接近精確解。綜合考慮計(jì)算精度和計(jì)算效率,下面的計(jì)算均取N=8。下面研究分別固有頻率與傾角和集中質(zhì)量的關(guān)系(過程中,出于貼近工程實(shí)際的考慮,任取水的流速U=35 m/s),臨界流速與傾角和集中質(zhì)量的關(guān)系,且如無特殊說明,所有參數(shù)均采用無量綱的形式。
下面計(jì)算當(dāng)γ=μ=0.01,且ξ0=0.4時(shí),傾角α變化時(shí)系統(tǒng)的前四階固有頻率隨傾角的變化關(guān)系,結(jié)果如表1所示。
表1 固有頻率隨傾角的變化
如表1所示,由于原始計(jì)算數(shù)據(jù)的原因,導(dǎo)致當(dāng)傾角α變化時(shí),各階固有頻率的變化并不明顯,但依然能夠發(fā)現(xiàn):當(dāng)管路向下傾斜(即α<π)時(shí),計(jì)算結(jié)果比向上傾斜(即α>π)時(shí)的大,意味著同等參數(shù)水平下,向上傾斜的管路的剛度更低。
(1)固有頻率與γ和μ的關(guān)系
取α=0,ξ0=0.2,且其余參數(shù)不變,通過計(jì)算,得到固有頻率隨γ和μ的變化曲線,結(jié)果如圖2所示。
圖2 固有頻率隨γ和μ的變化曲線
由圖2發(fā)現(xiàn),隨著γ和μ的增加,所有前四階固有頻率在最初階段均急速減小,而后減小趨勢放緩。
(2)固有頻率與ξ0的關(guān)系
取α=0,γ=μ=0.01,且其余參數(shù)不變,通過計(jì)算,得到固有頻率隨ξ0的變化曲線,結(jié)果如圖3所示。
圖3 固有頻率隨ξ0的變化曲線
由圖3發(fā)現(xiàn),隨ξ0的增加,固有頻率產(chǎn)生波動(dòng)形式的變化,此外,由于此處管路兩端的支承形式均為固定式,導(dǎo)致計(jì)算結(jié)果關(guān)于ξ0=0.5對稱。
(1)臨界流速與γ和μ的關(guān)系
取α=0,ξ0=0.2,其余參數(shù)不變,通過計(jì)算,得到在不同組γ和μ的前提下固有頻率隨流速的變化曲線,結(jié)果如圖4所示。
圖4 不同γ和μ下固有頻率隨流速的變化關(guān)系
根據(jù)圖4可以發(fā)現(xiàn):①隨著γ和μ的增加,在未發(fā)生任何形式的失穩(wěn)前,所有固有頻率均減小;②當(dāng)γ和μ繼續(xù)增大,前兩階模態(tài)的發(fā)散失穩(wěn)臨界流速始終保持不變,其中,一階對應(yīng)u1=6.284,二階為u2=8.987。
(2)臨界流速與ξ0的關(guān)系
取α=0,γ=μ=0.01,且其余參數(shù)不變,通過計(jì)算,得到在不同ξ0的前提下固有頻率隨流速的變化曲線,結(jié)果如圖5所示。
圖5 不同ξ0下固有頻率隨流速的變化關(guān)系
由圖5結(jié)合圖4(a),發(fā)現(xiàn)無論ξ0如何變化,一階模態(tài)發(fā)散臨界流速始終保持u1=6.284不變,但在此臨界值以后出現(xiàn)的失穩(wěn)形式有所差異,具體表現(xiàn)為:當(dāng)ξ0=0.1時(shí),二階模態(tài)首先出現(xiàn)發(fā)散,且臨界流速為u2=8.987,之后在u=9.084開始一階、二階模態(tài)合二為一;而當(dāng)ξ0=0.2時(shí),僅二階模態(tài)發(fā)生發(fā)散失穩(wěn),且與之對應(yīng)的臨界流速仍為u2=8.987,之后直到u=10,一階、二階模態(tài)始終處于發(fā)散狀態(tài);但是,當(dāng)ξ0=0.5時(shí),由于一階模態(tài)發(fā)散的提前結(jié)束,二階模態(tài)在未發(fā)生發(fā)散時(shí)便與一階模態(tài)耦合,發(fā)生顫振,與之對應(yīng)的流速的臨界值為uc=8.976,之后二者合二為一,然后在u=9.121再次發(fā)散。
(1)通過引入重力和傾角這兩個(gè)因素,建立了含集中質(zhì)量的斜置輸流直管的橫向振動(dòng)微分方程,利用基于新型形函數(shù)的Galerkin法對該微分方程進(jìn)行了離散,推導(dǎo)得到了計(jì)算管路固有頻率的特征方程。
(2)經(jīng)過實(shí)例計(jì)算,發(fā)現(xiàn):同等參數(shù)水平下,向上傾斜的管路的剛度更低;隨集中質(zhì)量的增加,在既定的流速下,固有頻率首先急速減小,然后趨勢放緩,但集中質(zhì)量的改變不會(huì)影響發(fā)散臨界流速的計(jì)算結(jié)果;隨集中質(zhì)量安置坐標(biāo)的增加,在既定的流速下,固有頻率波動(dòng)變化,且安置坐標(biāo)能夠引起臨界流速和失穩(wěn)形式的變化。
(3)研究可被借鑒用于研究具有其他支承形式、附加元素的管路的流固耦合振動(dòng)問題,為后續(xù)的振動(dòng)可靠性等問題的研究作了良好的理論及計(jì)算方法方面的鋪墊。