基于ANCF的旋轉超彈性厚壁REF振動特性分析

范 博, 王忠民

(1. 西安理工大學 機械與精密儀器工程學院,西安 710048; 2. 西安理工大學 土木建筑工程學院,西安 710048; 3. 黃河科技學院 建筑工程學院, 鄭州 450006)

彈性基環模型(ring on the elastic foundation,REF)由彈性基和柔性環兩部分組成,被廣泛應用于柔性齒輪[1]、子午線輪胎[2]和機械彈性輪[3]等旋轉柔性部件的動力學特性分析中。其彈性基一般由人工彈簧組和內壓兩部分組成,通過模態分析試驗[4-5],可以得到人工彈簧組的周向/徑向剛度。根據柔性環的不同假設,平面內REF模型可分為三種類型:基于Euler-Bernoulli彎曲梁的經典REF、基于Timoshenko彎曲梁的薄壁REF和基于平面應力問題的厚壁REF。

經典的REF最初是由Tielking[6]提出的。在隨后續的研究中,將經典REF的柔性環視為由不可擴展線彈性材料[7]、可擴展線彈性材料或可擴展超彈性材料[8]制成的Euler-Bernoulli曲梁。考慮到橫向剪切應變,Vu等[9]假設柔性環為Timoshenko曲梁,建立了一個薄壁REF模型,并分析了其線性振動特性。Lu等[10]考慮了Timoshenko梁理論之外的高階剪切修正因子,建立了高階薄壁REF。岳曉峰等[11]結合試驗數據,利用有限元軟件建立了超彈性薄壁REF,分析了子午線輪胎的振動特性。在上述REF模型中,柔性環都被視為環形梁,故忽略了柔性環徑向應變對于模型動力學特性的影響,因此這些模型只能用于分析薄壁結構。在Lu等[12-13]的后續研究中,考慮了柔性環的徑向應變,并開發了一個線彈性厚壁REF模型來探討彈性基剛度對柔性環穩態響應的影響。就作者所知,關于超彈性厚壁REF的動力學研究目前在國內外仍然是一個空白。

Shabana[14]提出的絕對節點坐標法(absolute nodal coordinate formulation,ANCF)框架不區分隨動坐標系,在固定坐標系中建立的動力學方程中質量矩陣為常數矩陣,且不含有離心力項和Coriolis力項。魏永[15]在不區分隨動坐標系的框架下,利用達朗貝爾原理推導了高速電鋸片的離心力項,在忽略Coriolis力條件下,分析了旋轉速度對于振動頻率的影響。Chen等[16]利用隨動坐標系和固定坐標系的轉化關系,推導了平面內柔性梁彎曲過程中廣義離心力和Coriolis力。Fan等[17]建立了具有隨動坐標系的ANCF框架,利用自適應曲梁單元建立了一種經典REF模型,探究了忽略Coriolis力條件下子午線輪胎的振動特性。

Mooney[18]通過物質相變理論和大量試驗,利用Cauchy-Green變形梯度張量的基本不變量建立了超彈性材料應變能密度Mooney-Rivlin模型。Yeoh模型[19]也采用相似的表達形式對超彈性材料的應變能密度進行定義。同樣階數下,由于不考慮Cauchy-Green變形梯度張量的第二不變量的作用,Yeoh模型要比Mooney-Rivlin模型在形式上要簡單。在文獻[20-22]中試圖在ANCF的框架下,結合上述兩種應變能函數模型對超彈性材料的動力學特性進行研究。但是由于這些模型中都含有Cauchy-Green變形梯度張量的基本不變量的一次項,這會導致在無應變的狀態下系統的仍具有較大的偽初始廣義彈性力。在其他文獻中沒有提出偽初始廣義彈性力存在的問題,更沒有對這個偽初始廣義彈性力進行補償。

在具有隨動坐標系的ANCF框架下,本文提出RAAE(rotating ANCF annular-sector element)單元,建立厚壁REF模型。基于Yoeh本構模型推導了厚壁REF模型的非線性運動微分方程,并基于振動平衡位置建立了該模型的線性化運動微分方程。通過與其他文獻結果進行對比,驗證了該模型的有效性,并通過一系列控制參數變量的計算預測了該模型的動力學特性。與其他研究不同的是,在本文中討論偽初始廣義彈性力對于振動平衡位置的影響,并且討論在Coriolis力作用下厚壁REF模型振動頻率隨旋轉變化的規律。

1 位置場描述

圖1 厚壁REF模型和RAAE單元

如圖1(b)所示,RAAE單元的變形過程可以被分解為2個運動的合成:一個是RAAE單元隨著隨動坐標系旋轉由絕對初始位置剛體運動到相對初始位置;另一個是RAAE單元在隨動坐標系內由相對初始位置發生柔性變形到相對變形位置。

在O-XY中,在隨動坐標系O-XY中,RAAE單元的相對位置可以用用無量綱坐標ξ和η的插值函數表示為

X=A0+A1ξ+A2ξ2+A3ξ3+A4η+A5η2+

A6η3+A7ξη+A8ξ2η+A9ξη2+A10ξ2η2+

A11ξ3η+A12ξη3+A13ξ3η2+A14ξ2η3+A15ξ3η3

Y=B0+B1ξ+B2ξ2+B3ξ3+B4η+B5η2+

B6η3+B7ξη+B8ξ2η+B9ξη2+B10ξ2η2+

B11ξ3η+B12ξη3+B13ξ3η2+B14ξ2η3+B15ξ3η3

根據雙三階Hermite插值函數,RAAE單元相對初始位置場函數R0和相對變形位置場函數R可以被表示為

(1)

R=S(ξ,η)qij(t)

(2)

RAAE單元的相對變形廣義單元坐標qij(t)可以寫為

qij=[(Nij)T(Ni(j+1))T(N(i+1)j)T(N(i+1)(j+1))T]T(3)

(4)

(6)

其中,

(7)

(8)

(9)

(10)

式中,φi=2π/n。

RAAE單元的形函數S(ξ,η)為2×32矩陣可以被表示為

S=[S11ajS21bjS12ajbjS22

S31ajS41bj+1S32ajbj+1S42

S13ajS23bjS14ajbjS24

S33ajS43bj+1S34ajbj+1S44]

(11)

在形函數S(ξ,η)中Sαβ(α,β為下標識符表示1,2,3,4)可以被寫為

Sαβ=Sα(ξ)Sβ(η)I2×2

(12)

(13)

設厚壁REF模型在隨動坐標系O-XY中初始和變形后全局相對位置的廣義坐標分別為q0和q(t),通過引入布爾矩陣Bij,可以將RAAE單元相對初始廣義坐標和相對變形后廣義坐標表示為

(14)

qij=Bijq(t)

(15)

(16)

(17)

其中,

(18)

(19)

(20)

(21)

(22)

(23)

式中,下標?=ξ或η為矩陣或向量對坐標ξ或η求一階導數。

根據式(23),絕對變形位置場函數對時間的一階導數和二階導數可以表示為

(24)

(25)

2 廣義彈性力和廣義剛度矩陣

2.1 超彈性厚壁環體的廣義超彈性力和廣義剛度矩陣

(26)

其中,

(27)

(28)

令

(29)

假設模型軸向線段在變形前后的長度和方向都不發生變化,那么模型在空間中Cauchy-Green應變張量C可以被寫為

(30)

(31)

(32)

將不可壓應變能函數理解為與體積變形無關的等積(偏量)應變能,再添加僅與體積比μ有關的體積(靜水)應變能部分。為了避免計算過程中的體積閉鎖現象,超彈性材料的Yeoh模型的應變能密度函數被表示為

(33)

式中:I1=tr(C);μ2=I3=det(C);α為體積應變能參數。

必須指出的是由于應變能密度函數中存在Cauchy-Green變形梯度張量的基本不變量的一次項,當超彈性材料在初始狀態下(此時,I1=3和μ=1),會存在

(34)

根據式(33)可得,超彈性環的彈性力的虛功可以被表示為

δWR(q)=δqTQR(q)

(35)

(36)

式中:h為環的橫向(垂直于平面方向)高度;Aj為單元初始面積,Aj=aj(bj+bj+1)/2;QR為超彈性環體的廣義彈性力列陣。

超彈性環體的廣義剛度矩陣可以被表示為

在沒有任何柔性變形發生的相對初始位置場,超彈性環體的初始廣義彈性力應當為0。根據式(34)和式(36)可以推導,在相對初始位置場,超彈性環體具有偽廣義彈性力,即QR(q0)≠0。在文獻[23]中都沒有對這個偽初始廣義彈性力進行討論和補償。

2.2 彈性基的廣義彈性力和廣義剛度矩陣

厚壁REF模型的彈性基由人工彈簧組和內壓兩部分組成。

根據胡克定律,環內壁單位面積內人工彈簧組的應變能面密度函數可以寫成

(38)

式中,u和v分別為環體內緣邊界上任意一點沿初始徑向和環向量的位移。它們可以根據幾何關系表示為

(39)

根據式(39),式(38)可以被改寫為

(40)

人工彈簧組彈性力的虛功可以被表示為

δWS(q)=δqTQS(q)

(41)

(42)

式中,QS為人工彈簧組的廣義彈性力列陣。

人工彈簧組的廣義剛度矩陣可以被表示為

(43)

假設旋轉超彈性厚壁REF在變形過程中內部壓力和溫度不會發生變化。在此理想假設下,內壓的應變能密度函數可以寫成

UP=P0ΔV

(44)

基于軸向線段不變假設,厚壁REF內腔的微元體積變化可以被表示為

(45)

根據式(45),式(44)可以被改寫為

(46)

內壓在變形過程中的虛功可以被表示為

δWP(q)=δqTQP(q)

(47)

(48)

式中:QP為內壓的廣義彈性力列矩陣;b1為最內側單元內緣的周向長度。

內壓的廣義剛度矩陣可以被表示為

(49)

3 非線性動力學微分方程

基于法線段不變假設,厚壁REF的慣性力做的虛功可以表示為

(50)

將式(21)和式(25)代入式(50),可得

(52)

式中,M為超彈性厚壁REF模型的廣義質量矩陣。

虛功原理為

δWT+δWR+δWS+δWP=δWW

(53)

式中,δWW為廣義外力的功。

將式(50),式(35),式(41)和式(47)代入式(52),并對偽初始廣義彈性力進行補償,經過運算可得旋轉超彈性厚壁REF模型的非線性運動微分方程為

(54)

QE(q)=QR(q)-QR(q0)+QS(q)+QP(q)

(55)

式中,單下劃線和雙下劃線項分別為廣義Coriolis力和離心力項。

若超彈性厚壁REF模型以勻速轉動(角速度Ω)時,式(54)可以被轉化為

(56)

令qS為平衡狀態下系統的廣義坐標,則

(57)

通過牛頓迭代法可以求得qS。

(58)

根據泰勒一階展開式,可得

(59)

將式(57)和(59)代入式(58),可得

(60)

設δqS=Aejωt代入式(60)得

(61)

式中,A為特征向量。式(61)有非零解的充要條件為

(62)

式中:ω=[ω(0,0)ω(1,0)ω(0,1)Fω(0,1)B…ω(N,M)Bω(N,M)F],下標N為節圓數,M為節徑數,F為前行波,B為后行波。需要指出的是,q為的全局相對變形廣義坐標,固有頻率中ω(N,M)F,ω(N,M)B為隨動坐標系O-XY中的前后行波頻率。

(63)

4 數值計算與分析

基于岳曉峰等的試驗數據,建立超彈性厚壁REF模型,相關參數如表1所示。

表1 模型參數

在模型中添加向心的集中載荷F=3 750 N。如圖2所示,當單元分布為20×1和30×2時,模型加載點的局部變形幾乎完全相同。這說明單元劃分為20×1時,超彈性厚壁REF模型具有較好的收斂性。

圖2 不同單元分布的受載模型的平衡位置

根據Yeoh模型的本構參數,可得等效線彈性厚壁REF模型的剪切模量G*=2C10,彈性模量E*=α以及泊松比ν*=0.5。如圖3所示,在相同載荷條件和單元離散方式下,超彈性和等效線彈性模型具有相似的平衡位置。這間接說明了本文模型的有效性。

圖3 超彈性模型和等效線彈性模型的平衡位置

如圖4所示,偽初始廣義彈性力不進行補償時,超彈性厚壁REF模型的平衡位置處,環壁厚度減小,并且每一個單元都出現了相似的畸變。為了保證所求平衡位置的準確性,偽初始廣義彈性力應當要進行補償。

圖4 無初始彈性力補償時受載模型的的平衡位置

在自由狀態下,部分試驗振型和本文模型振型如圖5所示,自由振動的部分頻率分析結果如表2所示。通過對比試驗結果和超彈性薄壁REF模型的結果,說明本文所建立的超彈性厚壁REF模型可以有效地分析子午線輪胎的振動特性。

表2 自由振動的部分頻率

圖5 自由狀態下振型

除此之外,利用超彈性厚壁REF模型分別計算得出了(0,0),(1,0)和(0,1)階頻率和振型,如表3所示。

表3 (0,0),(1,0)和(0,1)階自由振動頻率

分別假設模型的彈性基參數為kξ=2.5×10λPa,kη=2×10λ2Pa和P0=2.5×10λ3Pa。如圖6所示,在表1所示模型參數的基礎上,通過改變單一彈性基參數可以發現:當彈簧剛度超過2×104Pa時,模型的各階頻率開始隨其的增大而增大。當徑向彈簧剛度超過2.5×1016Pa,切向彈簧剛度超過2×1012Pa模型的各階頻率開始趨于收斂。內壓大于2×104Pa時,模型的各階頻率會隨著內壓的增大而增大。由于內壓過大,不符合實際情況,故沒有進一步的分析。

圖6 彈性基參數對于振動頻率的影響

如圖7所示,由于模型環內緣被彈性基所約束,造成模型結構剛度分布不均勻。隨著環壁厚度a的增加,超彈性模型的各階頻率先增加后減少。當環壁厚度為0.05 m時,相對于其他組別,模型的各階頻率最大。高階頻率受環壁厚度影響比低階頻率更加明顯。

圖7 環壁厚度對于振動頻率的影響

圖8 在中,角速度對于振動頻率的影響

5 結 論

本文基于絕對節點坐標法,利用平面內旋轉環扇形單元構建了厚壁REF模型,并應用此模型分析了子午線輪胎的面內振動特性。具體結論如下:

(1)采用Yeoh模型和近似不可壓縮超彈性應變能密度模型時,由于其含有Cauchy-Green變形梯度張量的基本不變量的一次項,在初始狀態下,會存在偽廣義彈性力,需要對其進行補償。否則,會嚴重影響系統平衡狀態的求解,導致計算結果的錯誤。

(2)模型各階振動頻率會隨著彈性基參數的增大呈現出3個變化階段,最初振動頻率緩慢增大,然后迅速增大,最后增大速度降低振動頻率趨于收斂。隨著環壁厚度的增大,模型各階振動頻率先增大,后減小。

(3)當振型周向振動波數為0時,其對應頻率不會隨著旋轉速度的變化而出現明顯的改變。當振型周向振動波數大于0時,隨著轉速的增大,對應頻率會分裂為2個行波頻率。算例中,子午線輪胎的第一階臨界角速度為285 rad/s。