王之千，毛保全，馮帥，楊雨迎，朱銳

(陸軍裝甲兵學(xué)院 兵器與控制系, 北京 100072)

0 引言

黏彈性膠體材料屬于非牛頓流體，是一種半液半固黏稠狀高分子聚合物，具有良好的彈性、可壓縮性及流動性等優(yōu)異的特點(diǎn)。將其裝入特殊的密閉空腔機(jī)械結(jié)構(gòu)中，起到減振、緩沖、吸能等效果，該裝置即為黏彈性膠體阻尼緩沖器。黏彈性膠體在阻尼緩沖器中的啟動流形式根據(jù)緩沖器機(jī)械結(jié)構(gòu)可分為間隙流和孔隙流兩種。

由于黏彈性流體獨(dú)特的非線性，國內(nèi)外眾多學(xué)者根據(jù)不同的研究背景和黏彈性材料，選用不同的數(shù)學(xué)模型和本構(gòu)關(guān)系研究阻尼緩沖器中黏彈性流體的流動特性。研究主要分為兩類：一類為經(jīng)驗(yàn)或半經(jīng)驗(yàn)公式，這類研究公式推導(dǎo)基于大量試驗(yàn)數(shù)據(jù)，所推出的經(jīng)驗(yàn)公式根據(jù)研究對象的不同而不同，使用限定條件較多，適用范圍較為局限。另一類為經(jīng)典模型及其衍生模型，其中以Maxwell模型最為典型。尤其隨著分?jǐn)?shù)階導(dǎo)數(shù)的提出，分?jǐn)?shù)Maxwell模型被證明可以較好地模擬黏彈性流體的流動過程[1-5]。Yang等[6]采用分?jǐn)?shù)Maxwell模型研究了黏彈性流體在無限長管道中恒定壓力梯度驅(qū)動下從靜止開始的單向啟動流動情況，求出了精確解。Yao[7]在逆行階梯Maxwell模型基礎(chǔ)上通過張量和不同分?jǐn)?shù)階導(dǎo)數(shù)定義將分?jǐn)?shù)階黏壺模型用于模擬大變形非線性黏彈性流體的流動。上述研究從數(shù)學(xué)角度很好地描述了黏彈性流體的流動特性，但值得注意的是，這些表達(dá)式模型參數(shù)往往通過數(shù)值方法得出，并不始終在物理和工程上保持合理性[8]。

圖1 彈簧壺元件原理圖Fig.1 Schematic diagram of a spring-pot as an element

黏彈性膠體材料在阻尼緩沖器中運(yùn)動時呈現(xiàn)非牛頓流體特性，具有較強(qiáng)的黏彈性，可采用分?jǐn)?shù)Maxwell模型進(jìn)行研究。但上述文獻(xiàn)中均求得模型的解析解，然而分?jǐn)?shù)階微積分方程求解析解難度較大、過程較為復(fù)雜，計算機(jī)模擬時運(yùn)算速度較慢；而數(shù)值解仿真模擬易于實(shí)現(xiàn)。并且采用以準(zhǔn)態(tài)特性為參數(shù)的分?jǐn)?shù)Maxwell模型模擬黏彈性流體在圓管內(nèi)和平板間加速流動的啟動流鮮有報道，同時在工程上用于研究阻尼緩沖器節(jié)流間隙和阻尼孔的分?jǐn)?shù)Maxwell啟動流模型也報道較少。

因此，本文針對黏彈性膠體在阻尼緩沖器阻尼孔和節(jié)流間隙中的啟動流問題，構(gòu)建以準(zhǔn)態(tài)特性為參數(shù)的分?jǐn)?shù)Maxwell模型，采用有限差分法求得模型的數(shù)值解，并設(shè)計啟動流測試試驗(yàn)系統(tǒng)，驗(yàn)證分?jǐn)?shù)Maxwell模型模擬黏彈性膠體阻尼緩沖器啟動流的準(zhǔn)確性。

1 分?jǐn)?shù)Maxwell啟動流模型構(gòu)建

1.1 分?jǐn)?shù)Maxwell模型

廣義黏彈性模型是將彈簧與黏壺串聯(lián)看作一個機(jī)構(gòu)元件，其本構(gòu)關(guān)系可寫為

(1)

采用Maxwell模型串聯(lián)兩個彈簧壺機(jī)構(gòu)元件(κ,α)和(ψ,β)(見圖2)，其中κ、ψ分別為不同準(zhǔn)態(tài)特性(Pa·sα、Pa·sβ)。由Maxwell模型關(guān)系可得

σa=σb=σ,

(2)

εa+εb=ε，

(3)

式中：σa和εa分別為彈簧壺a元件(κ,α)的應(yīng)力與應(yīng)變；σb和εb分別為彈簧壺b元件(ψ,β)的應(yīng)力與應(yīng)變；σ和ε分別為Maxwell模型的總應(yīng)力與總應(yīng)變。

圖2 分?jǐn)?shù)Maxwell模型Fig.2 Fractional Maxwell model

將(1)式代入(2)式、(3)式，得

(4)

即分?jǐn)?shù)Maxwell模型為

(5)

1.2 結(jié)構(gòu)及簡化

某型黏彈性膠體阻尼緩沖器結(jié)構(gòu)如圖3所示。阻尼緩沖器支座用于連接阻尼緩沖器和安裝機(jī)構(gòu)，起一定承力作用。將黏彈性膠體材料通過黏彈性膠體注入器注入阻尼緩沖器腔室中，根據(jù)需要調(diào)整黏彈性膠體材料的填充量，使阻尼緩沖器具有一定的預(yù)壓力，稱為初始預(yù)緊力。當(dāng)阻尼緩沖器受到外界壓力或沖擊大于初始預(yù)緊力時，活塞桿推動活塞壓縮阻尼緩沖器腔室內(nèi)的黏彈性膠體材料。膠體材料被迫流過節(jié)流間隙和活塞上的阻尼孔，產(chǎn)生黏滯力，阻礙活塞前行，這一過程中外力轉(zhuǎn)化為熱能和勢能，消耗和存儲。當(dāng)外力被撤消時，沖擊力小于初始預(yù)緊力，黏彈性膠體材料自行膨脹，釋放壓縮過程中存儲的勢能，將活塞推回到初始位置，而通過阻尼孔和節(jié)流間隙被擠入前腔體的膠體材料流回后腔體，等待下次沖擊或壓力。

圖3 某型黏彈性膠體阻尼緩沖器示意圖Fig.3 Schematic diagram of a viscoelastic elastomer shock absorber

為符合黏彈性膠體阻尼緩沖器啟動流時的實(shí)際工況，將黏彈性膠體在阻尼緩沖器阻尼孔和節(jié)流間隙中啟動流的運(yùn)動關(guān)系簡化至二維空間中，分別看作流體在圓管內(nèi)加速流動(見圖4(a)，即流體在兩平行板間以U=at的速度流動，而兩平行板始終保持不動)和兩平行板間一板不動、一板加速運(yùn)動這兩種形式(見圖4(b)，即下板以U=at的速度向右移動，上板始終保持不動，流體由于黏性在下板的帶動下向右流動)，圖中L為阻尼孔或節(jié)流間隙間距，U=at表示流速，a為加速度。

圖4 阻尼緩沖器結(jié)構(gòu)與黏彈性膠體的運(yùn)動關(guān)系Fig.4 Motion relationship between shock absorberstructure and viscoelastic elastomer

1.3 控制方程

流體運(yùn)動中總動量的時間變化率等于其自身體積力和相互作用表面力的總和。將加速度看作單位質(zhì)量流體動量隨時間的變化率，根據(jù)動量平衡率的牛頓運(yùn)動定律可知：

Ma=∑F,

(6)

式中：M為質(zhì)量；F為體積力和表面力總和。

與坐標(biāo)系無關(guān)的運(yùn)動方程矢量表達(dá)式為

(7)

式中：f為單位體積力；ρ為黏彈性流體密度；σ為應(yīng)力張量。

考慮阻尼緩沖器中黏彈性流體介質(zhì)運(yùn)動時黏性較大、總質(zhì)量較小，體積力遠(yuǎn)小于黏性力，因此將體積力忽略不計，則有

(8)

流速只在y軸方向有分布，則

(9)

式中：u(t,y)和σ(t,y)分別表示在t時刻空間位置y處的流速和應(yīng)力。

(10)

(11)

由Riemann-Liouville分?jǐn)?shù)階微積分公式可知：

(12)

(13)

Γ(z)是伽馬函數(shù)。

1.4 初始值和邊界條件

當(dāng)t≤0時，u(t)=0；當(dāng)t>0時，阻尼孔啟動流中心點(diǎn)的流速以加速度a逐漸增大，而靠近管壁的流速始終為0 m/s. 在節(jié)流間隙中，靠近運(yùn)動板的流體流速隨運(yùn)動板以加速度a逐漸增大，而靜止板附近的流速始終為0 m/s. 因此，初始條件和邊界條件為

1)阻尼孔：

(14)

2)節(jié)流間隙：

(15)

1.5 無量綱分析

采用Π定理分析分?jǐn)?shù)Maxwell模型的無量綱變量，

(16)

將(16)式的無量綱變量分別代入(10)式、(11)式中，得到分?jǐn)?shù)Maxwell啟動流無量綱方程(為簡便起見，省略無量綱標(biāo)記“*”)：

(17)

(18)

對應(yīng)的阻尼孔和節(jié)流間隙的無量綱初始和邊界條件為

(19)

(20)

1.6 數(shù)值求解

(21)

(22)

式中：o(Δt)為皮亞諾余項(xiàng)。

則

(23)

式中：t0=0；h0=1；hi=(i+1)1-h-i1-h，i=0,1,…,k.

同理，利用中心差分定理：

(24)

根據(jù)分?jǐn)?shù)階積分定義：

(25)

(26)

(27)

(28)

另外，牛頓流體模型根據(jù)(1)式可寫為

(29)

(30)

2 數(shù)值計算結(jié)果及分析

將模型數(shù)值解和邊界條件代入數(shù)學(xué)仿真軟件MATLAB中，為保證結(jié)果收斂，選取時間步長Δt=0.001、空間步長Δy=0.025. 通過數(shù)值計算，對比分?jǐn)?shù)Maxwell模型和牛頓流體模型的流動特性，選取基本模型參數(shù)為α=0.8、β=0.3、ζ=0.2.

2.1 阻尼孔數(shù)值模擬

圖5 阻尼孔流速分布Fig.5 Flow velocity distribution in orifice at different times

圖5為分?jǐn)?shù)Maxwell模型和牛頓流體模型隨時間增大時阻尼孔流速的變化情況。從圖5中可以看出，隨著時間的不斷增加，分?jǐn)?shù)Maxwell模型流速分布曲線由扁平狀逐漸發(fā)展為拋物線狀，符合實(shí)際黏彈性流體的流動特性，且流速分布曲線保持著良好的非線性特點(diǎn)，曲線曲率較牛頓流體大。這意味著分?jǐn)?shù)Maxwell模型模擬的流體始終保持著較好的黏彈性。

圖6為阻尼孔中心點(diǎn)的流體應(yīng)力變化情況。隨著時間的增加，分?jǐn)?shù)Maxwell模型和牛頓流體模型的應(yīng)力均呈線性增長。由于分?jǐn)?shù)Maxwell模型是通過廣義黏彈性模型推導(dǎo)衍生而得，而牛頓流體模型則是廣義黏彈性模型的一種特殊形式，因此二者應(yīng)力變化趨勢基本相同。但分?jǐn)?shù)Maxwell模型應(yīng)力小于牛頓流體模型。

取t=0.8，比較參數(shù)變化對分?jǐn)?shù)Maxwell模型流速的影響，如圖7～圖9所示。參數(shù)α和β均為分?jǐn)?shù)Maxwell模型的階數(shù)，二者的變化情況對流場的影響基本一致，隨著α和β的增加，分?jǐn)?shù)Maxwell模型流體曲率逐漸下降，但變化不大。而ζ相比α和β對曲線曲率影響較大。從圖9中可以看出，當(dāng)ζ=0.1時流速曲率最大，而隨著ζ的增加，流速曲率逐漸減小，并向牛頓流體發(fā)展。ζ與準(zhǔn)態(tài)特性有關(guān)，而準(zhǔn)態(tài)特性與松弛時間有關(guān)，松弛時間越長，材料越接近彈性體。這意味著ζ越大，黏度越小。

圖7 不同α?xí)r分?jǐn)?shù)Maxwell模型阻尼孔流速分布Fig.7 Flow velocity distribution of fractional Maxwell model in orifice for different α

圖8 不同β時分?jǐn)?shù)Maxwell模型阻尼孔流速分布Fig.8 Flow velocity distribution of fractional Maxwell model in orifice for different β

圖9 不同ζ時分?jǐn)?shù)Maxwell模型阻尼孔流速分布Fig.9 Flow velocity distribution of fractional Maxwell model in orifice for different ζ

2.2 節(jié)流間隙數(shù)值模擬

阻尼緩沖器節(jié)流間隙啟動流流速分布隨時間的變化情況如圖10所示。從圖10中可以看出，隨著時間的增加，牛頓流體模型流速分布逐漸向直線形狀發(fā)展，而分?jǐn)?shù)Maxwell模型始終保持著良好的非線性。進(jìn)一步證明了阻尼孔流速模擬的結(jié)果，分?jǐn)?shù)Maxwell模型非線性良好，具有較強(qiáng)的黏彈性。

圖10 節(jié)流間隙流速分布Fig.10 Flow velocity distributions in gap at different times

節(jié)流間隙啟動流靠近移動平板一側(cè)(y=0)應(yīng)力分布曲線如圖11所示。由圖11可見，隨著時間的增加，分?jǐn)?shù)Maxwell模型與牛頓流體模型應(yīng)力均呈線性增長，該結(jié)果與阻尼孔啟動流應(yīng)力分布增長趨勢一致。當(dāng)t=6時，分?jǐn)?shù)Maxwell模型和牛頓流體模型節(jié)流間隙啟動流應(yīng)力值分別為1.33和1.52，而阻尼孔啟動流應(yīng)力值分別為2.17和2.54；分?jǐn)?shù)Maxwell模型和牛頓流體模型阻尼孔啟動流的應(yīng)力分別為節(jié)流間隙啟動流應(yīng)力的1.63和1.67倍。原因在于，阻尼孔啟動流最大速度發(fā)生在平板間距中心點(diǎn)處(y=0.5)，而節(jié)流間隙啟動流最大速度出現(xiàn)在靠近移動平板一側(cè)處(y=0)，二者邊界條件不同。當(dāng)時間和空間步長相同時，阻尼孔啟動流計算所得速度數(shù)據(jù)是節(jié)流間隙啟動流數(shù)據(jù)數(shù)量的一半。應(yīng)力通過速度差值計算所得，而分?jǐn)?shù)Maxwell模型和牛頓流體模型應(yīng)力均呈線性增長，致使阻尼孔啟動流應(yīng)力值遠(yuǎn)比節(jié)流間隙啟動流大。雖然減小步長可以使流速數(shù)據(jù)數(shù)量保持一致，但步長的調(diào)整會導(dǎo)致差值系數(shù)發(fā)生變化，因此計算結(jié)果將保持不變。

圖11 節(jié)流間隙靠近移動平板一側(cè)應(yīng)力分布Fig.11 Stress distribution in gap near the side of the moving plate

3 試驗(yàn)驗(yàn)證

試驗(yàn)系統(tǒng)由伺服電機(jī)、電機(jī)控制器、玻璃管緩沖器、高速相機(jī)、圖像采集平臺和射燈光源組成，如圖12所示。為了方便捕捉黏彈性流體流動過程，采用塑料管作為活塞，玻璃管作為緩沖器腔壁，并按照緩沖器活塞和阻尼孔實(shí)際尺寸設(shè)計試驗(yàn)裝置，其中腔體內(nèi)壁直徑20 mm，活塞厚度15 mm，阻尼孔直徑2.5 mm. 為方便拍攝阻尼孔的流動過程，取消活塞與腔壁之間的間隙，并在塑料管活塞頭部安裝密封圈，防止運(yùn)動過程中有黏彈性流體從縫隙中擠入。

圖12 啟動流試驗(yàn)測試系統(tǒng)Fig.12 Test system for the start-up flow

采用日本Mitsubishi公司生產(chǎn)的交流伺服電機(jī)作為速度源推動活塞運(yùn)動，通過電機(jī)控制器調(diào)整轉(zhuǎn)速，使活塞桿始終保持恒定速度或加速度。利用美國Fastec公司生產(chǎn)的HiSpec 5 高速相機(jī)和配套的圖像采集平臺捕捉活塞啟動時黏彈性流體在阻尼孔中的流動情況。其中，拍攝過程中采用射燈作為輔助光源以提高相機(jī)識別度。

在玻璃管緩沖器腔體中注入主體成分為甲基硅氧烷的黏彈性膠體材料，其動力黏度為5×105mm2/s，密度為1.1 kg/m3，調(diào)節(jié)電機(jī)控制器使活塞桿速度從0 mm/s啟動，加速度始終保持1.0 mm/s2，在常溫下進(jìn)行試驗(yàn)。高速相機(jī)采集頻率為500幀，記錄下不同時刻黏彈性流體流動情況，如圖13所示，圖13中黑色部分為黏彈性膠體材料，白色部分為塑料管活塞。試驗(yàn)過程中輔助光源在玻璃管緩沖器左側(cè)，高速相機(jī)拍攝時左半部分容易出現(xiàn)反光點(diǎn)，導(dǎo)致試驗(yàn)結(jié)果左側(cè)有部分區(qū)域出現(xiàn)缺失，變?yōu)榘咨?。由于加速度較小，前0.6 s流體速度也較小，膠體變化不明顯。

圖13 高速相機(jī)記錄下阻尼孔啟動流運(yùn)動情況Fig.13 Test results of start-up flow in orifice

采用表1所列模型參數(shù)，代入分?jǐn)?shù)Maxwell模型中得到阻尼孔啟動流流速隨時間的變化結(jié)果。對速度進(jìn)行積分，獲得分?jǐn)?shù)Maxwell模型流體的液面高度分布曲線，如圖14所示。對比觀察試驗(yàn)和模型結(jié)果可知：阻尼孔中心位置流速最大，而靠近阻尼孔邊緣位置流速較小；隨著時間和活塞速度的增大，試驗(yàn)和數(shù)值結(jié)果的變化趨勢相同，液面形狀均從扁平狀逐漸發(fā)展成為拋物線狀，且阻尼孔啟動流過程始終保持著較強(qiáng)的黏彈性。

表1 分?jǐn)?shù)Maxwell模型參數(shù)

圖14 分?jǐn)?shù)Maxwell模型阻尼孔流體液面高度分布Fig.14 Fluid level distribution of fractional Maxwell model in orifice

其中，松弛模量Gt、松弛時間λ和分?jǐn)?shù)階指數(shù)α、β通過動態(tài)模量理論公式與試驗(yàn)數(shù)據(jù)擬合而得(見文獻(xiàn)[20])；而準(zhǔn)態(tài)特性參數(shù)κ和ψ[21]通過上述參數(shù)近似獲得：

(31)

以t=1.0 s時黏彈性膠體液面分布試驗(yàn)結(jié)果為例，提出最大液面高度H和最大曲率距離d兩個考量參數(shù)，如圖15所示。最大液面高度為黏彈性流體啟動時初始位置到當(dāng)前時刻液面最高點(diǎn)的垂直距離，是檢驗(yàn)同一時刻模型結(jié)果與試驗(yàn)結(jié)果的一致性。最大曲率距離是從液面最高點(diǎn)到液面最低點(diǎn)做輔助線，平行于輔助線的最大液面切線到輔助線之間的距離，主要用來考察流體的黏彈性。根據(jù)實(shí)際標(biāo)尺，測量最大液面高度和最大曲率距離，結(jié)果如表2所示，其變化趨勢如圖16和圖17所示。

圖15 黏彈性流體液面測量方法Fig.15 Viscoelastic fluid level measurement method

t/sH/mmd/mm試驗(yàn)結(jié)果分?jǐn)?shù)Maxwell模型結(jié)果牛頓流體模型結(jié)果試驗(yàn)結(jié)果分?jǐn)?shù)Maxwell模型結(jié)果牛頓流體模型結(jié)果0.60.37780.36000.36000.16090.15600.07820.80.66120.63460.63460.21020.21780.10961.00.93300.97780.97780.26370.26410.13321.21.39911.39421.39420.30880.29640.14651.41.81511.88041.88040.32500.30980.1549

圖16 最大液面高度結(jié)果對比Fig.16 Comparison of measured and simulated maximum fluid levels

圖17 最大曲率距離結(jié)果對比Fig.17 Comparison of maximum curvature distances

結(jié)合圖16、圖17及表2可知：牛頓流體模型隨著時間的增大，最大曲率距離有趨向于穩(wěn)態(tài)的趨勢，模型的黏彈性很弱；而分?jǐn)?shù)Maxwell模型結(jié)果與試驗(yàn)結(jié)果吻合較好，其中最大液面高度最大誤差為4.79%，最大曲率距離的最大誤差4.67%，如圖18所示。由此可見，分?jǐn)?shù)Maxwell模型可以很好地描述黏彈性膠體的黏彈性特點(diǎn)，準(zhǔn)確地模擬黏彈性膠體在阻尼孔中啟動流的實(shí)際流動過程。

圖18 分?jǐn)?shù)Maxwell模型的誤差Fig.18 Errors of fractional Maxwell model

4 結(jié)論

本文根據(jù)黏彈性膠體在阻尼緩沖器阻尼孔和節(jié)流間隙中的實(shí)際流動情況，將其簡化為黏彈性流體在圓管內(nèi)加速流動和兩平行板間一板不動、一板加速運(yùn)動兩種單向加速啟動流形式，根據(jù)實(shí)際工作環(huán)境設(shè)定初始和邊界條件，提出并構(gòu)建了含準(zhǔn)態(tài)特性參數(shù)用于研究黏彈性膠體阻尼緩沖器啟動流的分?jǐn)?shù)Maxwell模型，采用有限差分方法求得了模型的數(shù)值解。通過與牛頓流體模型進(jìn)行對比，分析了流速和應(yīng)力的變化情況及各參數(shù)對流速的影響。最后，設(shè)計了試驗(yàn)系統(tǒng)，捕捉了啟動流時黏彈性膠體材料在阻尼孔中的流動過程。通過定義、測量和對比分?jǐn)?shù)Maxwell模型模擬結(jié)果與試驗(yàn)結(jié)果的最大液面高度和最大曲率距離驗(yàn)證了模型的準(zhǔn)確性。得出以下結(jié)論：

1)黏彈性流體在圓管內(nèi)和兩平行板間加速流動時，分?jǐn)?shù)Maxwell模型流速分布始終保持著良好的黏彈性。

2)分?jǐn)?shù)Maxwell模型與牛頓流體模型應(yīng)力呈現(xiàn)線性增長趨勢。

3)隨著參數(shù)α、β和ζ的增加，分?jǐn)?shù)Maxwell模型流速分布曲線曲率逐漸下降，其中ζ對曲線曲率的影響比α和β大。

4)啟動流試驗(yàn)測試系統(tǒng)捕捉到的黏彈性膠體在阻尼孔中啟動流形狀從扁平狀逐漸發(fā)展成為拋物線狀，且中心點(diǎn)處速度最大，流動過程具有較強(qiáng)的黏彈性。該結(jié)果與分?jǐn)?shù)Maxwell模型模擬結(jié)果變化趨勢相同。

5)分?jǐn)?shù)Maxwell模型模擬結(jié)果的最大液面高度和最大曲率距離與試驗(yàn)結(jié)果較為吻合，其最大誤差分別為4.79%和4.67%，證明了模型的準(zhǔn)確性。

因此，分?jǐn)?shù)Maxwell模型可準(zhǔn)確地模擬黏彈性膠體在阻尼緩沖器阻尼孔和節(jié)流間隙中啟動流的流動情況，可用于指導(dǎo)黏彈性膠體阻尼緩沖器的設(shè)計。