基于線性化假設(shè)的流體流動(dòng)仿真

基于線性化假設(shè)的流體流動(dòng)仿真

傳統(tǒng)的管道流固耦合分析通常采用cordic（1904）方法。Skalak(1956)發(fā)展了Joukouski的分析,包含了泊松耦合,Walker和Philips推導(dǎo)了6方程模型,含有徑向平衡與連續(xù)方程,Wiggert等(1985)簡(jiǎn)化了上述模型,推導(dǎo)出了泊松耦合的4方程模型。Budny(1991)引進(jìn)了粘性摩擦項(xiàng),即頻率相關(guān)摩擦項(xiàng)的影響,進(jìn)一步考慮了結(jié)構(gòu)阻抗的影響,并給出了特征線解法。在管道的頻率特性的計(jì)算中,求解任意兩點(diǎn)之間的傳遞矩陣是非常重要的,D′Souza和Oldenburger(1964)給出了管道中任意兩點(diǎn)之間流體的傳遞關(guān)系(純流體之間),對(duì)管道網(wǎng)絡(luò)特性分析的發(fā)展起到了重要的推動(dòng)作用。Lesmez(1990)推導(dǎo)了流固耦合系統(tǒng)的模態(tài)分析模型及傳遞矩陣,但沒有考慮摩擦項(xiàng)的影響,而這種摩擦項(xiàng)對(duì)諧振峰幅值有較大的影響。1流體的流變學(xué)模型流體與管道的耦合振動(dòng)包括節(jié)點(diǎn)耦合與泊松耦合。首先應(yīng)建立包括這兩種耦合在內(nèi)的振動(dòng)的數(shù)學(xué)模型。作以下假設(shè):①管道為等截面、純彈性、均質(zhì)和各向同性的圓管;②忽略管道徑向變形引起的流體徑向運(yùn)動(dòng)及流體繞管軸的旋轉(zhuǎn)運(yùn)動(dòng);③忽略管道的結(jié)構(gòu)阻尼;④流體的摩擦是與頻率相關(guān)的摩擦;⑤流體的流速遠(yuǎn)小于流體的音速;⑥流體內(nèi)不出現(xiàn)氣穴或體積分離現(xiàn)象。在管道與流體微元中,坐標(biāo)選擇如圖1所示。(1)fzt1ap-euzz-155555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555555?V?t+1ρf?Ρ?z+τ0πDAfρf=0?Fz?z1Ap+τ0πDAp-ρp?Uz?t=0?Ρ?t-2νΚ*?Uz?z+Κ*?V?z=0?Fz?t1Ap-E?Uz?z-νRe?Ρ?t=0?V?t+1ρf?P?z+τ0πDAfρf=0?Fz?z1Ap+τ0πDAp?ρp?Uz?t=0?P?t?2νK??Uz?z+K??V?z=0?Fz?t1Ap?E?Uz?z?νRe?P?t=0式中:K*=K/[1+2RK(1-ν2)/(eE)];τ0=4ρfνfR。(V-Uz)+2ρfνfR∫i0?(V-Uz)?t(β)W(t-β)?dβ;W是無窮小的時(shí)間間隔τ的函數(shù)τ=νft/R2;(2)明確xp、fffffffbb.zbtbfbfb.z-ry的計(jì)算?Μy?t=EΙp?Ry?zΜy?z+Fx=(ρpΙp+ρfΙf)?Ry?t?Fx?z=(ρpAp+ρfAf)?Ux?t?Fx?t=ApG(?Ux?z-Ry)(3)fyp、fyz?Μx?t=EΙp?Rx?z?Μx?z-Fy=(ρpΙp+ρfΙf)?Rx?t?Fy?z=(ρpAp+ρfAf)?Uy?t?Fy?t=ApG(?Uy?z+Rx)(4)管道傳遞矩陣的確定?Μz?z=ρpJ?Rz?t;?Rz?z=1GJ?Μz?t以上各式中用到的符號(hào)集中說明如下:V為流體流速;P為壓強(qiáng);U為管道速度;F為管道內(nèi)力;R(帶上、下標(biāo))為管道轉(zhuǎn)動(dòng)角速度;M為管道所受力矩;A為面積;I為慣性矩;J為管道橫截面極慣性矩;R為管道半徑(不帶上、下標(biāo));D為管道內(nèi)直徑;e為管道壁厚;ν為管材泊松比;E為管材彈性系數(shù);K為流體體積彈性模量;K*為經(jīng)過修正的流體體積彈性模量;ρ為密度;G為剪切模量;νf為流體運(yùn)動(dòng)粘度;S=jω;下標(biāo):p為管道;f為流體;上標(biāo):x,y,z分別表示x,y,z軸向。對(duì)以上4組方程進(jìn)行拉氏變換并整理用矩陣形式表示為于是合在一起組成一個(gè)14階矩陣方程˙X=[˙Xfp˙Xxz˙Xyz˙Xtz]14×1=[Τfp0000Τxz0000Τyz0000Τtz]14×14[XfpXxzXyzXtz]=ΤX解這個(gè)微分方程可得:X=AX0。X0為初始條件;A為傳遞矩陣。其中:A=EeλIzE-1,λ為T矩陣的特征值;E是各λ所對(duì)應(yīng)的特征向量所組成的矩陣。為了確定截面上的未知狀態(tài)量,通常要與出口截面的已知狀態(tài)量聯(lián)立,故需要倒推求解,對(duì)X=AX0作以下變換:X0=A-1x=Ee-λIzE-1X=GX。把G定義為管道的傳遞函數(shù),由以上推導(dǎo)可知:G是管長(zhǎng)z和S的函數(shù)。以上各式中,A(S)為頻率相關(guān)項(xiàng)τ0的Laplace的變換形式的系數(shù):τ0(S)=L[τ0(t)]=ρfRSθ1(j√S/νfR)-2(V-Uz)=A(S)(V-Uz)A(S)=ρfRSθ1(j√S/νfR)-2=ρfRSθ1(jλ)-2=ρfRS2(α-1)α=jλJ0(jλ)2J1(jλ);λ=√S/νfRα=0.5λ+0.25+316λ+316λ2+63256λ3式中:J0,J1分別為0和1階的第1類貝塞爾函數(shù)。為了計(jì)算方便,在G矩陣的推導(dǎo)中可先對(duì)幾個(gè)子矩陣分別進(jìn)行推導(dǎo),最后再綜合得到管道傳遞矩陣。設(shè):Tfp對(duì)應(yīng)的特征值為λ1～λ4,相應(yīng)的特征向量組成的矩陣為Efp;Txz對(duì)應(yīng)的特征值為λ5～λ8,相應(yīng)的特征向量組成的矩陣為Exz;Tyz對(duì)應(yīng)的特征值為λ9～λ12,相應(yīng)的特征向量組成的矩陣為Eyz;Ttz對(duì)應(yīng)的特征值為λ13和λ14,相應(yīng)的特征向量組成的矩陣為Etz。于是可得E=[Efp0000Exz0000Eyz0000Etz]?E-1=[E-1fp0000E-1xz0000E-1yz0000E-1tz]e-λΙz=[e-λ1z00??00e-λ2z0??000e-λ3z0?0???????????000??0e-λ14z]14×14至此可以通過矩陣相乘求得直管的傳遞矩陣:G=Ee-λΙzE-12管系傳遞矩陣的修正在液壓管路系統(tǒng)中處理局部管路數(shù)學(xué)模型時(shí)多采用局部坐標(biāo)系,但由于復(fù)雜管路系統(tǒng)的表達(dá)要在一個(gè)統(tǒng)一坐標(biāo)系中進(jìn)行,這樣就需要坐標(biāo)的變換。采用合適的局部坐標(biāo)系對(duì)處理系統(tǒng)在節(jié)點(diǎn)處的變換關(guān)系是至關(guān)重要的,作者根據(jù)需要采用如下的方式來建立局部坐標(biāo)系:①以流體管路的傳輸方向?yàn)閦軸方向;②對(duì)于彎曲的管路,總以相鄰兩管路所構(gòu)成的平面的法向?yàn)閤方向。根據(jù)以上原則,如圖2所示空間管系中一段管路的求解過程為:(1)確定彎管處i點(diǎn)的x坐標(biāo),i點(diǎn)的x軸方向應(yīng)為管段i-1—i和i—i+1所在平面的法向方向,可以通過如下公式求得xi=Ri-1,i×Ri,i+1(2)從iL點(diǎn)到i-1點(diǎn)的坐標(biāo)變換,由于z軸不變,只需調(diào)整x,y軸即可,調(diào)整角α為α=arccosxi-1xi|xi-1||xi|則變換矩陣為(3)從iR到iL的坐標(biāo)變換,經(jīng)過一彎管,彎角為由于x軸不變,則坐標(biāo)變換只是針對(duì)y,z坐標(biāo)進(jìn)行,變換矩陣為這里需要特別注意的是,流體與固體在該點(diǎn)存在節(jié)點(diǎn)流固耦合問題,如下各量還有一增量存在:以上各式中的上標(biāo)R和L分別表示左邊和右邊。可以通過上式修正傳遞矩陣中對(duì)應(yīng)的各項(xiàng)元素,得到節(jié)點(diǎn)的傳遞矩陣。這樣沿著一定的方向依次將各管段、各節(jié)點(diǎn)和各元件的傳遞矩陣相乘,就可以得到整個(gè)管系的傳遞矩陣了。采用本文提到的兩次坐標(biāo)變換,只需要知道管路的彎曲角度和該段管路與其前后管段所組成平面的夾角,就可以得出所需矩陣,具有很強(qiáng)的通用性。3空間管系的密度圖3所示的是整個(gè)實(shí)驗(yàn)系統(tǒng),可調(diào)速的泵控馬達(dá)用于驅(qū)動(dòng)液壓泵。圖4～圖6所示的是3種類型的空間管系,管道壁厚為1mm,半徑為8mm,油液彈性模量為1500MPa,密度為872kg/m3,管材彈性模量為210GPa,密度是7800kg/m3,泊松比為0.3,管路彎曲角度均為90°。(1)項(xiàng)的影響時(shí)幅值不正確從圖7和圖8的仿真曲線中可以看出,未考慮頻率相關(guān)摩擦項(xiàng)的影響時(shí),幅頻特性曲線諧振峰的幅值遠(yuǎn)遠(yuǎn)超出實(shí)際情況,并且諧振峰幅值并不隨頻率增加而改變,這顯然也是不正確的?？紤]頻率相關(guān)摩擦影響的仿真曲線明顯更接近實(shí)際情況。(2)比較模擬結(jié)果和實(shí)驗(yàn)結(jié)果實(shí)驗(yàn)中選擇測(cè)量管道應(yīng)力σ來驗(yàn)證仿真結(jié)果。圖9～圖11分別為不同情況時(shí)的仿真與實(shí)驗(yàn)對(duì)照曲線。實(shí)線為仿真曲線,虛線為實(shí)驗(yàn)曲線。4幅值仿真結(jié)果分析由圖7～圖11的仿真結(jié)果和實(shí)驗(yàn)結(jié)果的對(duì)比可以得出如下結(jié)論:(1)摩擦項(xiàng)對(duì)振動(dòng)幅值有重要影響,如果不考慮摩擦項(xiàng)的影響將給諧振峰幅值帶來相當(dāng)