王亞軍，陳鼎銘，方紅榮，賀啟林，周浩洋

(1. 中國航天電子技術(shù)研究院，北京 100094；2. 北京宇航系統(tǒng)工程研究所，北京 100076)

0 引言

網(wǎng)套補償器是航天管路系統(tǒng)中廣泛采用的管配件，其主要功能是減振、補償位移。在工程實際中，鋼絲網(wǎng)套補償器往往承受包括溫度載荷、位移載荷、內(nèi)壓載荷以及沖擊載荷等不同形式的載荷，是管路結(jié)構(gòu)中容易破壞的部位[1]。為設(shè)計高性能、高可靠性的鋼絲網(wǎng)套補償器，以滿足日益提高的工程需求，必須對其力學(xué)性能開展深入、全面的研究。

補償器力學(xué)分析的難點在于復(fù)雜的微結(jié)構(gòu)特征：網(wǎng)套補償器由波紋管與鋼絲網(wǎng)套焊接組合而成，鋼絲網(wǎng)套又由多錠鋼絲編織而成，受載時網(wǎng)套和波紋管之間將發(fā)生強烈相互作用，使得網(wǎng)套補償器的剛度表現(xiàn)出強烈的非線性[2]。

考慮到補償器中波紋管的相關(guān)分析方法已經(jīng)比較成熟，鋼絲網(wǎng)套的分析方法是目前的研究重點。盛冬平等[3]、張程[4]、楊燕等[5]在補償器仿真時使用基于節(jié)點耦合的梁單元建立了網(wǎng)套模型，評估了網(wǎng)套補償器的靜動力學(xué)特性?；诹簡卧木W(wǎng)套建模方法對結(jié)構(gòu)進行了較多簡化，計算效率高，但難以保證準(zhǔn)確性，且不能計算軸向變形時剛度較低的階段。一些學(xué)者使用實體單元將網(wǎng)套的細(xì)觀編織結(jié)構(gòu)引入仿真模型，Rial等[6]、Huang等[7]建立了網(wǎng)套的實體模型，計算了補償器的軸向拉伸剛度。網(wǎng)套實體模型比梁模型更符合實際結(jié)構(gòu)，但對計算資源的要求極高，并且在軸向拉伸量較大時難以收斂。另有學(xué)者嘗試建立金屬絲網(wǎng)套的均質(zhì)化本構(gòu)，Hachemi等[8]在非線性正交異性層合板本構(gòu)中引入了編織角變化，以體現(xiàn)網(wǎng)套拉伸過程中的非線性行為，但在拉伸曲線后半段偏差較大。胡牧原等[9]進一步引入了基于試驗數(shù)據(jù)的編織角加速系數(shù)以增強剛度非線性。

本文以網(wǎng)套補償器這一復(fù)雜工程元件為研究對象，通過簡化鋼絲網(wǎng)套模型，使用仿真方法研究了網(wǎng)套的軸向剛度特性。建立了基于接觸關(guān)系的補償器有限元仿真模型，實現(xiàn)了對網(wǎng)套補償器軸向拉伸全過程的靜剛度分析。

1 鋼絲網(wǎng)套軸向剛度分析

鋼絲網(wǎng)套結(jié)構(gòu)是補償器軸向剛度特性的非線性來源，也是分析難點。本文建立了網(wǎng)套的簡化分析模型——螺旋梁模型，通過有限元仿真方法計算了幾何參數(shù)對其軸向剛度特性的影響，研究了邊界條件對螺旋梁軸向剛度的影響。

1.1 網(wǎng)套簡化分析模型——螺旋梁模型

1.1.1 螺旋梁模型簡介

網(wǎng)套軸向變形時由于軸向?qū)ΨQ性和反對稱性，其中鋼絲的響應(yīng)在統(tǒng)計意義上是相同的，因此在軸向變形計算時，可以以單根鋼絲為分析對象。網(wǎng)套的細(xì)觀編織結(jié)構(gòu)使得鋼絲的路徑比較復(fù)雜，難以用簡單的數(shù)學(xué)公式描述。若忽略鋼絲的卷曲，則鋼絲的路徑可簡化為螺旋線(見圖1)，使分析難度大大降低。

圖1 螺旋梁模型示意圖Fig.1 Spiral beam model of wire cloth

簡化之后，鋼絲的路徑幾何由中心線的螺旋角α0、金屬絲螺旋直徑D0、網(wǎng)套軸向長度l0這3個幾何參數(shù)確定。對于由N錠、每錠m根鋼絲編織而成的網(wǎng)套，其軸向剛度Kz可由單根鋼絲軸向剛度Ks疊加得到

Kz=mNKs

(1)

使用空間梁單元建立螺旋梁仿真模型，對其進行軸向(圖1中z方向)變形響應(yīng)分析。作為分析算例的網(wǎng)套規(guī)格見表1。

表1 網(wǎng)套參數(shù)表Tab.1 Mesh sleeve parameters

1.1.2 螺旋梁模型可靠性驗證與誤差修正

采用文獻[10]中管狀編織結(jié)構(gòu)幾何建模方法建立表1中鋼絲網(wǎng)套的編織單胞仿真模型(見圖2)，分別對內(nèi)部無約束工況及內(nèi)部具有剛性芯工況下螺旋梁模型和單胞模型進行仿真分析，以驗證螺旋梁模型在模擬網(wǎng)套單獨拉伸及與內(nèi)部芯軸接觸時與單胞模型的誤差。

圖2 鋼絲網(wǎng)套編織單胞仿真模型Fig.2 Simulation model of steel wire mesh sleeve unit cell

用兩種仿真模型得到的兩種工況下網(wǎng)套軸向拉伸的力-位移曲線見圖3。

(a)無內(nèi)部約束

由圖3及表2可知，網(wǎng)套單獨拉伸時螺旋梁模型得到的軸向剛度結(jié)果和編織單胞差別較??；而在內(nèi)部有剛性芯時，螺旋梁得到的剛度要顯著大于編織單胞，此時有必要對螺旋梁進行修正，以使其軸向剛度與編織單胞相等。本文采用將螺旋角從α0修正至γ0的方式使其軸向剛度與編織單胞等效，按照此方式在內(nèi)部有芯軸時，表1中網(wǎng)套的編織角由48°增大至52.857°。

表2 鋼絲網(wǎng)套軸向剛度仿真結(jié)果對比Tab.2 Comparisons of axial stiffness simulation results of steel wire mesh sleeve

1.2 幾何參數(shù)對軸向剛度的影響

1.2.1 軸向長度l0

螺旋角α0、網(wǎng)套直徑D0以及鋼絲直徑d與表1中網(wǎng)套相同，建立不同軸向長度l0的螺旋梁有限元仿真模型。工程實際中的網(wǎng)套邊緣通常使用邊條焊接固定，在軸向變形時兩端徑向、環(huán)向以及轉(zhuǎn)動自由度都被固定，只存在軸向位移。因此仿真時螺旋梁端點邊界位移和轉(zhuǎn)角條件如下

Ur|A=Ur|B=0
Uθ|A=Uθ|B=0
Uz|A=0,Uz|B=Wz
URr|A=URr|B=0
URθ|A=URθ|B=0
URz|A=URz|B=0

(2)

仿真獲得的力-應(yīng)變曲線見圖4。曲線圖表明，不同長度的螺旋梁在軸向拉伸時力-應(yīng)變曲線差別十分明顯。相同軸向應(yīng)變的情況下，螺旋梁的初始軸向長度越長軸向力越小，即其剛度與軸向長度l0之間并非線性關(guān)系。

圖4 不同軸向長度下螺旋梁軸向拉伸力-應(yīng)變仿真曲線Fig.4 Axial tensile stress-strain simulation curves of helical beams with different axial lengths

對于特定長度的螺旋梁，其拉伸曲線在彈性范圍內(nèi)隨著軸向應(yīng)變增大，曲線斜率也隨之升高。這由幾何非線性導(dǎo)致，拉伸過程中梁的大變形(見圖5)影響了剛度。軸向長度較短的螺旋梁(如l0=7.5 mm時的曲線)在較大軸向力范圍內(nèi)剛度變化較小，軸向長度較長的螺旋梁(如l0=27.5 mm和l0=30 mm時的曲線)在較大的應(yīng)變范圍內(nèi)剛度變化也較小。而軸向長度適中的螺旋梁具有最顯著的幾何非線性。

圖5 l0=30 mm螺旋梁拉伸變形位移云圖(結(jié)果經(jīng)柱坐標(biāo)系下z軸旋轉(zhuǎn)陣列展示)Fig.5 Spiral beam tensile deformation displacement nephogram with l0=30 mm(the results are displayed by z-axis rotation array in cylindrical coordinate system)

從仿真結(jié)果中提取軸向應(yīng)變?yōu)?.001時的軸向反力，擬合初始的小變形范圍內(nèi)的軸向剛度(見圖6)。由圖6可知，螺旋梁處于較小軸向長度的范圍時，剛度對軸向長度較為敏感，此時剛度隨初始長度增加而急劇減?。宦菪旱某跏驾S向長度較大時，軸向剛度極低，并且剛度隨長度增加而減小的趨勢也放緩。

圖6 螺旋梁初始軸向拉伸剛度與軸向長度的關(guān)系Fig.6 Relationship between initial axial tensile stiffness and axial length of spiral beam

通過提取螺旋梁變形時的應(yīng)變能仿真值可以解釋造成上述現(xiàn)象的原因。在螺旋梁中總的應(yīng)變能中由拉伸變形貢獻的應(yīng)變能可表示為

(3)

其中，L為螺旋梁的總弧長，εN為沿梁的中心線切向應(yīng)變，σN為沿梁的中心線切向應(yīng)力。

由彎曲、扭轉(zhuǎn)和剪切變形貢獻的應(yīng)變能則為

EM=Etotal-EN

(4)

對于Euler-Bernoulli梁，計算時忽略了剪切變形，因此只存在拉伸應(yīng)變能和彎扭應(yīng)變能。從仿真結(jié)果提取的各軸向長度下，螺旋梁內(nèi)拉伸應(yīng)變能與總應(yīng)變能的比值見圖6中的紅色虛線，可見拉伸應(yīng)變能占總應(yīng)變能比例越高，螺旋梁的軸向剛度也相應(yīng)越大。這表明螺旋梁剛度隨軸向長度的非線性變化主要是內(nèi)部變形模式的轉(zhuǎn)變導(dǎo)致的。螺旋鋼絲軸向長度較大時拉伸變形很小，主要依靠彎扭實現(xiàn)軸向伸長，此時剛度很低。而軸向長度很小時，軸向拉伸時鋼絲內(nèi)主要發(fā)生拉伸變形，剛度要高得多。

1.2.2 螺旋角α0

軸向長度l0、網(wǎng)套直徑D0以及鋼絲直徑d與表1中網(wǎng)套相同，建立不同螺旋角α0的螺旋梁仿真分析模型。仿真獲得的力-應(yīng)變曲線如圖7所示。曲線圖表明，螺旋角越小，螺旋梁軸向剛度越大，同時拉伸曲線越接近線性。

圖7 不同螺旋角下螺旋梁軸向拉伸力-應(yīng)變仿真曲線Fig.7 Simulation curves of axial tensile force strain of spiral beam under different spiral angles

初始的小變形范圍內(nèi)(軸向應(yīng)變0.001)的軸向剛度仿真值如圖8所示。

圖8 螺旋梁初始軸向拉伸剛度與螺旋角的關(guān)系Fig.8 Relationship between initial axial tensile stiffness and helix angle of spiral beam

由圖8可知，螺旋梁軸向剛度對初始螺旋角較為敏感，二者之間為非線性變化關(guān)系，初始螺旋角越大剛度越低。通過對應(yīng)變能的計算可以發(fā)現(xiàn)，軸向剛度的下降與拉伸應(yīng)變能比例的下降趨勢基本一致，因此螺旋梁剛度隨初始編織角非線性變化的主要原因同樣是變形模式的轉(zhuǎn)變。

1.2.3 網(wǎng)套直徑D0

軸向長度l0、螺旋角α0以及鋼絲直徑d與表1中網(wǎng)套相同，建立不同網(wǎng)套直徑D0的螺旋梁仿真分析模型。仿真獲得的力-應(yīng)變曲線如圖9所示。曲線圖表明，相對于其他參數(shù)，網(wǎng)套直徑對螺旋梁軸向剛度影響較小。網(wǎng)套直徑越小，螺旋梁軸向剛度越小。

圖9 不同直徑下螺旋梁軸向拉伸力-應(yīng)變仿真曲線Fig.9 Axial tensile force strain simulation curves of spiral beam with different diameters

初始的小變形范圍內(nèi)(軸向應(yīng)變0.001)的軸向剛度仿真值如圖10所示。

圖10 螺旋梁初始軸向拉伸剛度與網(wǎng)套直徑的關(guān)系Fig.10 Relationship between initial axial tensile stiffness of spiral beam and diameter of mesh sleeve

由圖10可知，螺旋梁軸向剛度隨著網(wǎng)套初始直徑D0線性變化，剛度隨初始直徑增大而緩慢升高。對應(yīng)變能的計算表明，網(wǎng)套直徑對拉伸應(yīng)變能比例影響較小，與對剛度的影響趨勢一致。

1.3 邊界條件對軸向剛度的影響

1.3.1 固定邊界與約束徑向位移的循環(huán)邊界

假設(shè)一種理想情形：拉伸軸向長度極大的補償器時，由于網(wǎng)套的徑向收縮效應(yīng)(見圖6)，其中的網(wǎng)套將與波紋管接觸。波紋管的表面形狀復(fù)雜，如果波紋管的剛性很大，則此時網(wǎng)套相當(dāng)于受到周期性的徑向約束，約束周期如圖11中AB段所示。

圖11 網(wǎng)套受波紋管周期性徑向約束示意圖Fig.11 Schematic diagram of periodic radial constraint of the mesh sleeve by the bellows

將上述問題進一步理想化為一根無限長度的螺旋梁受到周期性的徑向位移約束。在一個約束周期內(nèi)，梁的兩個截面應(yīng)滿足循環(huán)邊界條件，以使截面處的內(nèi)力和變形滿足連續(xù)性。在計算螺旋梁z方向拉伸響應(yīng)時，循環(huán)邊界條件要求端點滿足如下約束方程

Ur|A-Ur|B=0
Uθ|A-Uθ|B=0
Uz|A-Uz|B=Wz
URr|A-URr|B=0
URθ|A-URθ|B=0
URz|A-URz|B=0

(5)

其中下標(biāo)A，B分別代表梁的兩個端點(見圖1)。

由于端點的轉(zhuǎn)動自由度未被完全約束，循環(huán)邊界下的螺旋梁所受約束要弱于固定邊界，因此結(jié)構(gòu)剛度應(yīng)該不大于固定邊界下的螺旋梁。這里對不同軸向長度l0的螺旋梁仿真模型施加約束徑向位移的循環(huán)邊界，并將仿真結(jié)果與之前固定邊界下的仿真結(jié)果進行對比。不同邊界下仿真獲得的力-應(yīng)變曲線如圖12所示。

圖12 不同邊界條件下螺旋梁軸向拉伸仿真力-應(yīng)變曲線Fig.12 Axial tensile simulation force strain curves of spiral beam under different boundary conditions

仿真結(jié)果顯示，對于不同軸向長度的網(wǎng)套，循環(huán)邊界和固定邊界下的差別都極小。在初始變形階段，兩種邊界的力-應(yīng)變曲線基本重合，只有在大變形時，曲線才略微有差別。這表明在計算螺旋梁軸向剛度時這兩種邊界條件可以互換，而不會對精度產(chǎn)生影響。

1.3.2 徑向位移不完全約束

現(xiàn)實中與網(wǎng)套連接或發(fā)生接觸的不可能是完全的剛體，彈性體之間的接觸將沿接觸面法向發(fā)生一定量的位移?？紤]到這一點，在螺旋梁拉伸仿真時在兩端同步施加不同比例的徑向收縮位移，以模擬螺旋梁與不同彈性體作用后對軸向剛度的影響。

建立螺旋梁仿真分析模型，軸向長度l0、螺旋角α0、網(wǎng)套直徑D0以及鋼絲直徑d都與表1中網(wǎng)套相同。計算時軸向位移載荷相同，設(shè)置不同的徑向位移(其值為負(fù))。在不同的徑向收縮量下的軸向拉伸仿真結(jié)果對比如圖13所示。

(a)仿真力-應(yīng)變曲線

由圖13可見，徑向收縮量越大軸向剛度越小，相對于幾何參數(shù)的影響效果，端點徑向收縮量對螺旋梁整體軸向剛度徑向位移的影響較小。

2 補償器軸向剛度分析

補償器的軸向剛度特性表現(xiàn)出高度非線性，本研究基于螺旋梁模型，提出了基于接觸關(guān)系的子網(wǎng)套剛度分析方法。在子網(wǎng)套分析方法的基礎(chǔ)上，建立了適用于不同軸向拉伸階段的螺旋梁-波紋管模型，同時給出了低剛度階段和高剛度階段的補償器剛度理論估計方法。

2.1 補償器試驗軸向剛度特性

王亞軍等[2]對網(wǎng)套補償器的軸向剛度特性進行了實驗研究，文獻中試驗件1在無內(nèi)壓狀態(tài)下的軸向拉伸曲線如圖14所示。

圖14 補償器軸向拉伸試驗曲線Fig.14 Compensator axial curves of tensile test

圖14中曲線表現(xiàn)出顯著的非線性特征，可劃分為3個階段。

1)低剛度階段：在拉伸初始階段剛度近似為常值，力-位移曲線近似為線性；

2)過渡階段：拉伸至一定位移后力-位移曲線表現(xiàn)出強烈的非線性，剛度開始加速上升；

3)高剛度階段：隨著補償器拉伸位移繼續(xù)增加，剛度增速減緩，力-位移曲線又有回歸線性的趨勢。

將文獻[2]中無內(nèi)壓和充壓狀態(tài)下的補償器拉伸、壓縮力-位移曲線繪于同一圖(見圖15)中進行對比。由圖15可知，對于未充壓狀態(tài)下的補償器壓縮，其曲線斜率(剛度)與拉伸的低剛度階段是相同的。充內(nèi)壓之后補償器相當(dāng)于初始施加了一個壓力推力，其結(jié)果是力-位移曲線發(fā)生了平移，但形狀基本不變。因此，只要獲得補償器未充壓狀態(tài)下的力-位移曲線，便可以進一步推導(dǎo)得出其余軸向工況下的曲線。

圖15 無內(nèi)壓和充壓狀態(tài)下補償器軸向力-位移曲線Fig.15 Axial force-displacement curves of compensator with/without inner pressure

2.2 基于接觸關(guān)系的子網(wǎng)套剛度分析方法

補償器軸向拉伸時表現(xiàn)出的非線性剛度特性由鋼絲大變形這一幾何非線性過程以及網(wǎng)套和波紋管接觸這一邊界非線性過程共同導(dǎo)致。結(jié)合1.3節(jié)中對螺旋梁力學(xué)特性的討論，本文提出了基于接觸關(guān)系的子網(wǎng)套剛度分析方法，可以綜合考慮兩種非線性因素。

補償器松弛狀態(tài)下網(wǎng)套半徑大于波紋管，而受到軸向拉伸時二者將發(fā)生接觸。當(dāng)網(wǎng)套與波紋管發(fā)生接觸之后，截取兩波峰之間的網(wǎng)套(見圖16)，被截取的網(wǎng)套兩端徑向位移被波紋管限制，顯然該段網(wǎng)套受到的約束程度要比循環(huán)邊界(圖16綠色區(qū)域)更強，而比固定邊界(圖16紅色區(qū)域)更弱。

圖16 補償器中的網(wǎng)套分段及對應(yīng)邊界條件Fig.16 Mesh sleeve segmentation and corresponding boundary conditions in compensator

1.3.1節(jié)已經(jīng)證明了在計算螺旋梁模型的軸向剛度時這兩種邊界條件差別可以忽略，則介于二者之間的真實邊界可以用這兩種邊界來替代。這樣就可以將發(fā)生接觸的網(wǎng)套分段，并將該段網(wǎng)套的邊界替代之后，將計算對象縮小為一小段網(wǎng)套，即子網(wǎng)套。

在拉伸初始階段網(wǎng)套和波紋管之間沒有接觸，補償器剛度為二者的軸向剛度疊加，因此試驗曲線中該階段表現(xiàn)為線性。當(dāng)拉伸到一定量之后，網(wǎng)套和波紋管之間的接觸關(guān)系使得軸向剛度呈現(xiàn)出強烈的非線性，剛度急劇上升。網(wǎng)套變形和接觸區(qū)域的變化如圖17所示。

圖17 網(wǎng)套變形與接觸示意圖Fig.17 Deformation and contact diagrams of mesh sleeve

接觸發(fā)生之后，網(wǎng)套中間部分(圖17中紅色環(huán)線)徑向位移受到波紋管約束，此時可以把原網(wǎng)套視作兩個長度只有原來一半的子網(wǎng)套；隨著軸向拉伸量的增大，波紋管接觸區(qū)域由中間向兩邊擴展，網(wǎng)套被進一步拆分為更多軸向長度與波距相同的很短的網(wǎng)套(中間段子網(wǎng)套)，以及末端兩個稍長的網(wǎng)套(邊緣段子網(wǎng)套)。

設(shè)波紋管與網(wǎng)套相接觸的波數(shù)為n，顯然整體網(wǎng)套的軸向力Fz、位移Wz與各子網(wǎng)套之間有如下關(guān)系

Fd=Fc=Fz
2Wd+nWc=Wz

(6)

其中，下標(biāo)d表示邊緣段子網(wǎng)套對應(yīng)的物理量，下標(biāo)c表示中間段子網(wǎng)套對應(yīng)的物理量。

則網(wǎng)套整體軸向剛度Kz可表示為

(7)

根據(jù)式(7)，在拉伸的初始階段，n=0，Kz=Kd，需要計算的只有單個邊緣段子網(wǎng)套(網(wǎng)套本身)。當(dāng)網(wǎng)套和波紋管接觸之后，n逐漸變大，根據(jù)1.2.1節(jié)中軸向長度對剛度的影響，Kd隨軸向長度不斷變短而急劇上升，此時剛度Kc基本不變，網(wǎng)套整體剛度急劇非線性上升，這一過程中邊界非線性是剛度變化的主要因素。當(dāng)接觸區(qū)域不再變化，即網(wǎng)套與波紋管穩(wěn)定接觸之后，n和邊緣段子網(wǎng)套初始長度不再改變，但剛度Kd和剛度Kc受幾何變形影響而增大，網(wǎng)套整體剛度仍在緩慢上升，此時影響網(wǎng)套剛度的只有幾何非線性因素。

2.3 低剛度階段的補償器軸向剛度分析方法

在低剛度階段網(wǎng)套與波紋管無接觸，補償器剛度等于二者的線性疊加。網(wǎng)套對應(yīng)螺旋梁的軸向剛度以及波紋管的軸向剛度可通過有限元仿真獲得，使用的模型與過渡階段相同，下文詳述。

2.4 過渡階段的補償器軸向剛度分析方法

過渡階段網(wǎng)套和波紋管之間的相互作用力水平較弱，此階段網(wǎng)套剛度的劇烈增加由網(wǎng)套接觸區(qū)域變化導(dǎo)致，分析重點為接觸關(guān)系變化過程。本文通過建立2/N波紋管-螺旋梁復(fù)合有限元仿真模型，結(jié)合插值獲得了該階段補償器的軸向力-位移曲線。

建立波紋管2/N(N仍表示網(wǎng)套鋼絲錠數(shù))模型(見圖18)及單根鋼絲螺旋梁模型，并在波紋管環(huán)向邊界處施加循環(huán)邊界條件。模型中螺旋梁使用實體單元，波紋管使用殼單元。螺旋梁和波紋管之間設(shè)置無摩擦“硬”接觸，對波紋管和螺旋梁兩端面施加相同的軸向位移。

圖18 2/N波紋管-螺旋梁復(fù)合模型Fig.18 2/N bellows with spiral beam composite model

根據(jù)1.1.2節(jié)，補償器在網(wǎng)套與波紋管接觸前采用未修正螺旋角的螺旋梁模型計算誤差較小，在網(wǎng)套與波紋管完全接觸之后，應(yīng)該使用修正過螺旋角的螺旋梁。分別建立螺旋角為α0和γ0的2/N波紋管-螺旋梁復(fù)合模型(后文簡稱α0模型和γ0模型)。從有限元結(jié)果中提取螺旋梁和2/N波紋管的位移及反力結(jié)果，按如下關(guān)系進行疊加可得補償器整體軸向剛度K

(8)

其中，Kb為2/N波紋管的軸向剛度；式中右側(cè)第一項為鋼絲網(wǎng)套的整體軸向剛度，第二項為整個波紋管的軸向剛度。

在接觸過程中過渡階段的拉伸曲線，可以基于α0模型和γ0模型的仿真數(shù)據(jù)進行插值獲得，插值方法如下。

設(shè)α0模型得到的仿真位移-力曲線函數(shù)關(guān)系為x1=f1(y)，γ0模型得到的仿真位移-力曲線函數(shù)關(guān)系為x2=f2(y)。α0仿真模型鋼絲和波紋管開始接觸時為插值起點，此時軸向拉力值為a；γ0仿真模型接觸波數(shù)不再變動為插值終點，此時軸向拉力值為b。構(gòu)造插值函數(shù)x=p(y)，插值函數(shù)需滿足

p(a)=f1(a)
p(b)=f2(b)

(9)

按以上條件可構(gòu)造插值公式

(10)

最終由式(10)獲得的y=p-1(x)即為該階段補償器拉伸響應(yīng)的力-位移關(guān)系。

需要指出的是，α0模型顯然可以得到接觸發(fā)生之前補償器的拉伸曲線，即低剛度階段曲線，只需將接觸發(fā)生之前的螺旋梁和2/N波紋管的仿真結(jié)果代入式(8)即可。

2.5 高剛度階段的補償器軸向剛度分析方法

補償器拉伸到一定位移后，網(wǎng)套和波紋管接觸關(guān)系達到穩(wěn)定，力-位移曲線進入高剛度階段，可以按2.2節(jié)中的網(wǎng)套分段方法分別對邊緣段子網(wǎng)套和中間段子網(wǎng)套進行仿真分析。

2.5.1 邊緣段子網(wǎng)套

使用三結(jié)點二次空間梁單元建立邊緣段子網(wǎng)套的螺旋梁模型。邊緣段子網(wǎng)套兩端邊界更接近固定邊界(見圖16)，故在仿真計算時對螺旋梁兩端施加固定邊界條件。

2.5.2 中間段子網(wǎng)套

進入高剛度階段后，中間段子網(wǎng)套與波紋管相互作用較為強烈?？紤]到波紋管徑向彈性變形將對網(wǎng)套剛度產(chǎn)生影響，本研究在仿真建模時，中間段子網(wǎng)套與邊緣段子網(wǎng)套使用了不同的建模方法。對應(yīng)的螺旋梁兩端施加循環(huán)邊界，并且建立了波紋管模型，以體現(xiàn)二者的相互作用。

中間段子網(wǎng)套對應(yīng)的分析對象為兩端施加循環(huán)邊界條件且與兩波峰接觸的螺旋梁，如圖19左圖所示。左圖仿真模型的接觸區(qū)域被人為地分割，增加了接觸狀態(tài)判斷的計算量，使得有限元計算收斂困難。為避免上述問題，在實際計算時波紋管采用圖19右圖所示單波模型。該模型與所要分析的中間段網(wǎng)套以及對應(yīng)的波紋管波谷段在幾何上不同，但在施加循環(huán)邊界條件后，經(jīng)過陣列和旋轉(zhuǎn)操作，都可變換為長度無限的帶網(wǎng)套的補償器，二者描述的是同一整體力學(xué)模型的一個周期，故其軸向力學(xué)響應(yīng)是相同的。

圖19 中間段子網(wǎng)套和波紋管模型的等效變換Fig.19 Equivalent transformation of the model of the intermediate sub mesh sleeve and bellows

對于圖19所示的仿真模型，可根據(jù)環(huán)向周期性進一步將模型縮小為原來的2/N(見圖20)，只需對波紋管左右邊界施加環(huán)向周期性條件即可。圖20為最終計算中間段子網(wǎng)套軸向響應(yīng)使用的仿真模型。

圖20 2/N單波-單錠螺旋梁復(fù)合模型Fig.20 2/N single wave with single spindle spiral beam composite model

復(fù)合模型中鋼絲使用C3D8R實體單元，波紋管使用S4R殼單元。鋼絲和波紋管之間設(shè)置無摩擦“硬”接觸，對波紋管兩端施加軸向位移載荷，同時通過螺旋梁兩端面約束的參考點施加相同的軸向位移。

3 算例與結(jié)果分析

3.1 補償器參數(shù)

本文以文獻[2]中的試驗件1為算例，在下文分析補償器的軸向剛度。其中波紋管的規(guī)格如表3所示。

表3 分析所用補償器中波紋管規(guī)格Tab.3 The specifications of bellows in the compensator analyzed

網(wǎng)套參數(shù)與表1相同。

3.2 低剛度階段及過渡階段計算結(jié)果

根據(jù)1.1.2節(jié)，對于本文研究的起始編織角為48°鋼絲網(wǎng)套，對應(yīng)的螺旋梁在接觸波紋管之前α0取48°，與波紋管穩(wěn)定接觸之后γ0取52.9°。

接觸發(fā)生之前α0模型得到的低剛度階段仿真力-位移曲線與試驗曲線進行對比，見圖21。在該階段仿真值和試驗值一致性較好。

圖21 低剛度階段仿真與試驗力-位移曲線對比Fig.21 Comparison of force-displacement curves between simulation and test at low stiffness stage

圖22給出了仿真模型的接觸應(yīng)力云圖，α0模型計算得到鋼絲和波紋管開始接觸時，軸向位移Wz為0.69 mm；γ0模型計算得到的接觸波數(shù)不再變動時，軸向位移Wz為2.32 mm。由此可確定過渡階段的起點和終點。

圖22 2/N波紋管-螺旋梁復(fù)合模型接觸應(yīng)力云圖Fig.22 Contact stress nephogram of 2/N bellows with spiral beam composite model

按式(10)擬合的插值曲線與試驗曲線對比如圖23所示，二者的一致性較好。由力-位移曲線可知，當(dāng)補償器開始拉伸時，試驗值與螺旋角為α0模型的仿真值更為接近。當(dāng)網(wǎng)套與波紋管開始接觸之后，α0模型誤差逐漸增大，試驗曲線向螺旋角γ0模型的仿真曲線靠近。本文提出的基于螺旋梁仿真值的插值方法可有效模擬補償器過渡階段的力-位移曲線。

圖23 基于仿真結(jié)果的過渡階段插值曲線Fig.23 Interpolation curves in transition phase based on simulation results

3.3 高剛度階段計算結(jié)果

按照子網(wǎng)套分析方法，穩(wěn)定接觸之后需將網(wǎng)套分解為如表4所列的子網(wǎng)套。

表4 子網(wǎng)套參數(shù)表Tab.4 Parameters of sub mesh sleeve

將中間段網(wǎng)套及邊緣段網(wǎng)套仿真結(jié)果按式(8)的關(guān)系整合后，與波紋管力-位移曲線疊加，獲得高剛度階段的補償器仿真力-位移曲線。截取補償器試驗力-位移曲線的高剛度部分(按2.3節(jié)中的計算結(jié)果為438 N之后的曲線)，其與仿真值對比如圖24所示。

圖24 高剛度階段仿真與試驗力-位移曲線對比Fig.24 Comparison of force-displacement curves between simulation and test at low stiffness stage

由圖24可見，仿真模型得到的高剛度階段力-位移曲線呈現(xiàn)弱非線性，與試驗結(jié)果一致性較好。

表5給出了由試驗曲線中高剛度階段近似線性部分的軸向剛度、仿真曲線近似線性部分的軸向剛度對比。仿真值與試驗值偏差較小，由本文提出的高剛度階段仿真計算方法是可靠的。

表5 高剛度階段補償器剛度對比Tab.5 Stiffness comparison of compensator in high stiffness stage

3.4 補償器軸向計算結(jié)果總結(jié)

將各階段仿真曲線相結(jié)合，結(jié)果如圖25所示。

圖25 補償器軸向拉伸全過程力-位移曲線仿真值與試驗值Fig.25 Comparison of compensator force-displacement curves between simulation and test in the whole process of axial tension

由圖25可知，仿真得到的軸向拉伸全過程力-位移曲線與試驗曲線一致性較好，能夠復(fù)現(xiàn)補償器的非線性剛度特性，證明了本文提出的補償器軸向響應(yīng)分析方法的可靠性。

4 結(jié)論

本文提出了基于有限元仿真的網(wǎng)套補償器軸向拉伸全過程的剛度計算方法，主要結(jié)論如下。

1)根據(jù)網(wǎng)套軸向變形特點將分析對象簡化為螺旋梁，使用有限元仿真分析了幾何參數(shù)對軸向剛度的影響，結(jié)果表明螺旋梁軸向剛度對軸向長度、螺旋角較為敏感，軸向剛度隨參數(shù)非線性變化的主要原因為鋼絲內(nèi)變形模式的轉(zhuǎn)變。

2)分析了不同邊界條件對螺旋梁軸向剛度的影響，對于軸向拉伸工況，固定邊界和循環(huán)邊界計算結(jié)果基本一致，同時徑向約束由剛性變?yōu)閺椥詴r軸向剛度將下降。

3)建立了基于接觸關(guān)系的子網(wǎng)套剛度分析方法，根據(jù)接觸狀態(tài)將網(wǎng)套分解為兩波峰之間的中間段子網(wǎng)套和兩端的邊緣段子網(wǎng)套，解釋了補償器拉伸過程中低剛度階段、過渡階段和高剛度階段的軸向剛度變化。

4)使用有限元仿真方法分別建立了適用于低剛度階段及過渡階段的2/N波紋管-螺旋梁復(fù)合模型、適用于高剛度階段的2/N單波-單錠螺旋梁復(fù)合模型，實現(xiàn)了對補償器低剛度階段、過渡階段及高剛度階段的全過程軸向響應(yīng)計算。算例表明仿真獲得的力-位移曲線與試驗曲線一致性較好。