Kevlar繩索非對稱遲滯模型及參數識別

易 琳,王 班,郭吉豐

（浙江大學電氣工程學院,浙江杭州310027）

Kevlar繩索非對稱遲滯模型及參數識別

易 琳,王 班,郭吉豐

（浙江大學電氣工程學院,浙江杭州310027）

針對Kevlar繩索準靜態加載下出現的非對稱遲滯效應,建立能夠靈活控制遲滯環形狀的推廣Bouc-Wen模型以描述高度非對稱的遲滯現象.結合數學模型與試驗結果特點提出適用于該模型的分步參數識別方法,通過提出的識別方法得到試樣模型參數,采用數值計算方法求解模型得到在循環載荷作用下的理論遲滯曲線.結果表明,理論計算值與試驗值較吻合,證明了模型的正確性與識別方法的有效性.

遲滯效應；Bouc-Wen模型；Kevlar繩索；參數識別

合成纖維是一種具有高比模量、高比強度和低伸長率的工程材料,由其制作的繩索結構在航空航天、機械、電力和電子工業中得到越來越廣泛的應用.其中以Kevlar和Twaron為代表的芳香族聚胺纖維［1］有著優異的機械性能、熱性能和化學性能,非常適合于航空航天領域應用,如抗空間碎片碰撞結構、降落傘繩、空間繩系系統等［2-3］.

在對采用芳綸繩索的空間繩系衛星系統進行動力學分析及控制時,考慮衛星剛體動力學模型的同時,還須考慮連接衛星繩索的力學模型.以往的研究多以彈簧阻尼模型來描述繩索的力學行為［4］,該類模型相對簡單,有利于降低控制器設計的復雜性和提高計算效率,但在描述繩索的一些非線性特性時與實際相差較大.Cheng等［5］在對單根Kevlar KM2纖維的機械特性研究中,采用實驗方法測得了本構關系在縱向、橫向以及扭轉方向的5個獨立參數,對由多根纖維加捻編織而成的繩索特性研究有一定的參考意義.Shim等［6］研究由Twaron纖維制作的盔甲面料動態力學性能,發現材料在高應變速率與低應變速率作用下相比強度和彈性模量均顯著增大,采用黏彈性開爾文模型解釋了這一現象.Chailleux等［7］提出芳綸紗線的非線性黏彈性和黏塑性拉伸行為的數學模型,指出紗線在小形變范圍內具有非線性應力-應變關系,描述了恒應力下紗線的拉伸蠕變規律和最大應力相關的塑性變形正向積累規律.遲滯效應是指在載荷作用下材料的受力和形變之間表現出具有記憶特征的遲滯時變關系,Bouc于1967年提出一種光滑遲滯模型,后經Wen［8］改進得到能夠產生系列不同遲滯曲線的模型,用以描述大量存在于機械、土木、地震和材料等工程中的遲滯現象.筆者在對Kevlar繩索的準靜態受力—形變關系進行實驗研究時發現,繩索不僅表現出蠕變、非線性剛度和阻尼、塑性變形等非線性現象,而且在周期往復加減載時表現出明顯的遲滯效應.這給繩索動力學模型的建立帶來了新的挑戰,然而關于Kevlar繩索的遲滯非線性特性的研究鮮有文獻報道.

Kevlar繩索遲滯效應產生的機理較復雜,不僅在微觀上與材料本身的分子結構有關,而且在宏觀上與纖維加捻捻度及纏繞并線方式有關,因此很難從遲滯效應產生機理上來建立遲滯特性模型.本文對Kevlar繩索的軸向恢復力-形變關系進行試驗研究,得到靜態加載方式下的恢復力-形變遲滯曲線,對測得曲線的觀察分析發現,其為非對稱遲滯曲線.基于原始Bouc-Wen對稱模型的滯環形狀分析,提出一種推廣Bouc-Wen模型描述非對稱遲滯現象.采用適合該模型的分步識別方法識別出模型中的各參數值,根據識別得到的參數代入模型進行數值計算,與試驗曲線進行對比,以驗證建立模型的準確性及參數識別方法的有效性.

如無特別說明,本文所述繩索形變及恢復力均沿繩索軸向,在準靜態加載條件下僅研究縱向形變與恢復力的關系.

1 Kevlar繩索非對稱遲滯模型

圖1 Kevlar繩索試樣實物Fig.1 Specimen of Kevlar tether

Kevlar材質繩索如圖1所示,繩索由3股紗線并線纏繞而成,每股紗線由若干根纖維單絲加捻而成,試樣直徑為1.5 mm,總長度為18.1 m,斷裂拉力約為1 300 N.未經加載的Kevlar繩索在循環加載后會產生一定的殘余形變,該形變在卸載后無法恢復,且殘余形變與加載歷史中的最大負載呈正相關,即最大負載越大殘余形變量越大,Chailleux等［7,9］對Kevlar材料的這種特性進行研究.其成因在于芳族高分子間形成的分子鏈在拉力作用下發生旋轉,導致材料發生塑性變形,此外,繩索結構內部存在間隙,經過拉力拉伸后結構內部間隙變小,且該變化不可逆.該現象類似于填充橡膠材料的馬林斯效應［10］.為了研究Kevlar繩索的遲滯特性,須消除塑性變形的影響,故在試驗之前對試樣施加合適的預載荷,保持足夠長時間直至繩索的形變量不再增加.試驗通過試樣一端固定,另一端懸掛系列不同質量重物的方式測量周期載荷作用下繩索形變量與恢復力的關系.測量的難點主要在于如何非接觸地精確測量微小形變,避免測量設備的摩擦力作用干擾對繩索施力.試驗裝置采用SensoPart公司生產的FT 80 RLA-500-S1L8激光位移傳感器測量繩索伸長量,傳感器量程為500 mm,分辨率為0.5 mm.通過RS485通訊將采集的位移數據傳輸至計算機存儲和處理,試驗裝置原理如圖2所示.圖中,x為形變,F為恢復力.試樣在進行遲滯特性測量試驗前經過持續12 h的85 N軸向預加載,所有試驗均在恒定室溫（20℃）的條件下進行.

圖2 Kevlar繩索遲滯特性試驗裝置原理圖及試驗結果Fig.2 Experiment set and result on hysteresis effect of Kevlar tether

1.1 Bouc-Wen模遲滯形狀控制

試驗結果表明,Kevlar繩索在準靜態加載下的軸向受力和形變關系具有較明顯的遲滯現象,表明繩索的形變不僅與受力有關,還與受力歷史過程有關,是典型的非線性時變系統,不能用簡單的非線性彈簧模型描述.Bouc-Wen模型具有通用性及便于數學處理等優點,因而被廣泛地應用于描述非線性遲滯系統,物理模型如圖3所示.

該模型能夠從現象上解釋和反應局部記憶遲滯效應,如下所示：

圖3 恢復力物理模型Fig.3 Physical model of restoring force

式中：t為時間,F（t）為繩索恢復力,Fa（t）為彈性恢復力,x（t）為形變量；z（t）為遲滯恢復力,是為反映遲滯特性引入的虛擬狀態變量,無實際物理意義,也無法直接測量；A和n分別為控制遲滯環大小和光滑程度的常數,ψ為關于x、和z的非線性函數,控制遲滯環形狀因素。通過嘗試不同非線性函數,對圖2的繩索恢復力和形變關系曲線進行擬合發現,曲線與三次多項式函數曲線較接近,令

z（t）不僅與當前形變量x（t）有關,而且和加載歷史過程有關；Fa（t）僅與x（t）有關,與加載歷史過程無關.

對于原始的Bouc-Wen模型,有

式中：γ和β為與Bouc-Wen遲滯環形狀有關的常數.

式（1b）等號兩邊同時乘以d t／d x,可得

從式（3）可以看出,x-z相平面內的遲滯環斜率d z／d x由模型參數A、n及形狀控制函數ψ決定,通常函數ψ的值僅與x、和z的符號有關.根據時間的物理意義有d t＞0,所以與d x的符號相同.根據式（3）可知,特別是當n≠1時,通常難以給出z關于x的解析形式,因此通過數值計算分析是有效的.

圖4 形狀控制函數值Fig.4 Shape control function

1.2 推廣Bouc-Wen模型

繩索只能承受軸向拉力,不能承受軸向壓力,相應的形變只能是拉伸形變（x＞0）.觀察圖2可以發現,遲滯曲線不存在與形變量x相關的幾何中心,從物理意義上來說遲滯特性不以某形變值為分隔.考慮形狀控制函數與x無關,對Bouc-Wen形狀控制函數改進來得到在4個相位內對應不同數值的形狀控制函數是亟待解決的.按照該思路,提出4自由度推廣Bouc-Wen模型：

式中：β1、β2、β3和β4為常數.該形狀控制函數包含4個相位,依據˙x和z的不同符號組合可以有4個不同的取值,圖4（b）表示模型在一個加載周期內的非對稱遲滯環.當ψ1＝ψ3且ψ2＝ψ4時,改進的非對稱遲滯模型退化為Bouc-Wen對稱模型,因此可以認為該模型是Bouc-Wen模型的一種推廣模型.圖4（b）中4個相位內ψ的取值與˙x和z的符號組合對應的關系如表1所示.

表1 推廣Bouc-Wen模型ψ取值Tab.1 Values ofψfor generalized Bouc-Wen model

表1列出的ψi與βi（i＝1,2,3,4）之間的線性關系可由下式表示：

式（5）中的系數矩陣滿秩,故有

利用式（5）、（6）的矩陣方程,可以通過非線性最優化方法,擬合試驗數據得到相應的ψi,最優目標為所有試驗數據與理論計算值的差值平方和最小,再由式（6）計算得到βi.直接根據式（4）和試驗數據可以識別出使試驗數據和理論計算結果最接近的一組模型參數βi,而無需考慮ψi.前一種識別方法的優點是參數更易于識別,特別是在分段識別各相位內的ψ時.

2 非對稱遲滯模型參數識別

2.1 最小二乘識別法

假定試驗測得的離散數據為（xi,Fi）（i＝1,…, m）,由式（1）和（4）組合形成了Kevlar繩索的非對稱遲滯恢復力模型.該模型包括9個獨立參數,形成以下向量：

將式（1a）、（1c）代入式（3）得到模型參數識別的誤差函數：

為了得到一組使理論模型計算結果最接近試驗數據的模型參數,可以采用最小二乘法進行模型參數識別.將數據點（xi,Fi）（i＝1,…,m）代入式（8）,得到對應的殘量：

目標函數為所有數據點的殘差平方和,即

式中：r（y）為殘量函數,r（y）＝［r1（y）,r2（y）,…, rm（y）］T.使式（10）最小的向量y是識別得到的參數,可以利用高斯-牛頓方法［11］求解式（10）的非線性最小二乘問題.

2.2 參數分步識別方法

吳善躍等［12］在研究鋼絲繩減震器的非對稱遲滯模型參數識別問題時提到,在非線性最小二乘問題求解過程中會出現矩陣病態導致識別效果不理想,提出一種具有更好識別效果的分解識別方法.本文研究的對象恢復力模型有所不同,且繩索結構有特殊性,如只能拉伸不可壓縮,故須探索適合Kevlar繩索非對稱遲滯模型參數識別效果更理想的方法.采用降低待識別參數向量y的維數的方法,根據參數類別結合參數在模型中所起的作用逐一識別,上一識別環節的結果作為下一環節的已知參數,直至得到參數向量中所有元素的識別結果.

2.2.1 彈性恢復力模型參數識別 從遲滯恢復力z的物理意義并根據式（3）可知,遲滯恢復力在x-z相平面有界,令d z／d x＝0,有

這表現在圖2（b）中遲滯回線隨x的增大逐漸接近上界zmax＝,反之逐漸接近下界zmin＝.從圖2中觀察到不同幅值的循環加載曲線在加載區域和減載區域各有一共同的漸進線,這2根漸近線的部分曲線段如AB和CD段所示,為遲滯恢復力z接近上、下界時形成.由式（1a）、（11）可知：

設（xABi,FABi）（i＝1,…,r）和（xCDi,FCDi）（i＝1,…,s）分別為AB和CD段內的試驗數據,將數據點寫成向量形式：xr＝［xAB1,xAB2,…,xABr］,Fr＝［FAB1,FAB2,…,FABr］,xs＝［xCD1,xCD2,…, xCDs］,Fs＝［FCD1,FCD2,…,FCDs］,代入式（12）有

式（13）完成了對試驗數據彈性恢復力部分與遲滯恢復力部分的解耦,經過處理后的數據可以直接用來識別彈性恢復力模型式（1c）.

將圖2中相同橫坐標x對應的加載部分FDAB與減載部分FDCB作減法運算得到的最大差值為18 N,可令zmax＝9,zmin＝－9.采用線性最小二乘識別方法：

彈性恢復力模型識別結果如圖5所示.

2.2.2 遲滯恢復力模型參數識別 由于遲滯恢復力無法直接測量得到,由式（1a）可知,將圖5中的恢復力實測值與彈性恢復力計算值相減,可得遲滯恢復力,即zi＝Fi－＋a1xi）,計算結果如圖6所示,將得到的遲滯恢復力計算結果作為遲滯恢復力測量值.

將z＝0代入式（3）,有d z／d x＝A,即根據圖6中曲線在z＝0處的斜率可以求得模型參數A.對于圖6中3次不同循環加載幅值下的遲滯恢復力曲線,分別求得相位4、1交越處與相位2、3交越處的斜率,再對求得的斜率求平均得到A＝380.由式（11）可知：

圖5 彈性恢復力模型參數識別結果Fig.5 Identification result of elastic restoring force

圖6 遲滯恢復力測量值Fig.6 Measured values of hysteretic restoring force

參數向量y中剩下3個待識別的元素：n、ψ2和ψ4.由圖6的數據（xi,zi）（（i＝1,…,m）和式（15）,結合y中已識別得到的參數代入式（8）～（10）,導出關于3個待識別元素的非線性最小二乘問題,采用高斯-牛頓迭代方法求解得到參數識別結果：n＝2.93, ψ2＝－6.543,ψ4＝－0.113.將n代入式（15）計算得到ψ1＝0.552,ψ3＝0.552.由式（6）計算可得：β1＝－1.388,β2＝1.608,β3＝－1.608,β4＝1.940.

2.2.3 模型參數識別過程總結 識別過程如圖7所示.首先要保證試驗數據覆蓋足夠大的形變范圍,以盡可能使遲滯恢復力達到“飽和”并且上下“飽和”區域對應的橫坐標x上有重合.觀察恢復力模型可以發現,彈性恢復力項只與形變x有關,與加載歷史過程無關,遲滯恢復力項在試驗結果中隨x表現出單調變化且存在上界和下界.利用這一性質,可以計算出遲滯恢復力的最大、最小值zmax和zmin.根據遲滯曲線的漸近線段部分數據可以識別彈性恢復力系數a1、a2和a3.將遲滯恢復力項分離出來,遲滯恢復力為系統遲滯現象的根源,根據d z／d x在z＝0處的值識別參數A,而ψ1、ψ3可以用參數n表示.然后采用非線性最小二乘識別方法,得到參數n、ψ2和ψ4.最后根據參數βi與ψi（i＝1,2,3,4）之間的線性關系計算出βi的值,完成對模型中所有參數的識別.

圖7 推廣Bouc-Wen模型參數識別過程Fig.7 Process of identification for generalized Bouc-Wen model

3 模型及其參數識別結果校驗

將參數識別結果代入非對稱遲滯模型式（1）和（4）,可得

針對圖2所示的測試系統,需要考慮重物自身的動力學模型：

式中：m為重物質量,試驗時通過改變重物質量來實現加載或減載；g為重力加速度,取9.8 m／s2.

以重物質量m作為系統輸入,重物重力周期變化且幅值在每2個周期減小一次,以模擬實際加載過程,如圖8所示.對式（16）、（17）進行數值計算,模擬在循環加載條件下F與x間的變化關系,仿真計算結果與試驗結果對比如圖9所示.圖9表明,采用提出的非對稱遲滯模型能夠反映試樣的遲滯特性,利用該模型來描述Kevlar繩索在靜態加載下的恢復力-形變關系是適合的；模型的仿真計算值與實測值比較吻合,證明提出的分步參數識別方法能夠根據循環加載試驗結果有效地識別出模型參數.

圖8 非對稱遲滯仿真模型輸入Fig.8 Simulation input of asymmetric hysteresis model

圖9 非對稱遲滯模型仿真計算結果與試驗結果Fig.9 Simulation result of asymmetric hysteresis model compared with experimental data

4 結 語

Kevlar繩索在工程應用中表現出強非線性,特別是大尺度的強能量耗散.本文針對其準靜態加載情況下軸向恢復力與形變關系表現出的非對稱遲滯現象,基于對工程領域應用廣泛的Bouc-Wen模型形狀控制分析,提出一種推廣的Bouc-Wen模型.該模型包括4個獨立的形狀控制參數,能夠分別控制遲滯環4個不同相位內的形狀,靈活地描述各種不同形狀的非對稱遲滯現象.提出適合Kevlar繩索特性的參數分步識別方法.將參數識別的結果代入數學模型,得到相應的計算值,與實測值相比較吻合.證明采用該模型描述Kevlar繩索的非對稱遲滯特性是合適的,提出的識別方法能夠有效地識別該模型的各參數,可為類似的芳綸材料編織繩帶結構遲滯特性研究提供參考.

［1］YANG H H.Kevlar aramid fiber［M］.New York：Wiley,1993.

［2］SPENCER D A,BLANCHARD R C,BRAUN R D,et al.Mars Pathfinder entry,descent,and landing reconstruction［J］.Journal of Spacecraft and Rockets,1999, 36（3）：357- 366.

［3］劉濱濤,賈光輝,黃海.Kevlar層合板超高速撞擊數值建模及參數識別［J］.宇航學報,2011,32（02）：261- 266.

LIU Bin-tao,JIA Guang-hui,HUANG Hai.Numerical modeling and parameter identification for Kevlar laminate under the condition of hypervelocity impact［J］.Journal of Astronautics,2011,32（02）：261- 266.

［4］SEBASTIAN M S,UNNIKRISHNAN K C,NARAYANAN S.Viscoelastic properties of Kevlar-29 fabric tape strength member［J］.Mechanics of Materials, 2008,40（11）：949- 960.

［5］CHENG M,CHEN W N,WEERASOORIYA T.Mechanical properties of Kevlar（R）KM2 single fiber［J］.Journal of Engineering Materials and Technology：Transactions of the ASME,2005,127（2）：197- 203.

［6］SHIM V P W,LIM C T,FOO K J.Dynamic mechanical properties of fabric armour［J］.International Journal of Impact Engineering,2001,25（1）：1- 15.

［7］CHAILLEUX E,DAVIES P.Modelling the non-linear viscoelastic and viscoplastic behaviour of aramid fibre yarns［J］.Mechanics of Time-Dependent Materials, 2003,7（3／4）：291- 303.

［8］WEN Y K.Method for random vabration of hysteretic systems［J］.Journal of the Engineering Mechanics Division-ASCE,1976,102（2）：249- 263.

［9］CHENG M,CHEN W N.Modeling transverse behavior of Kevlar（R）KM2 single fibers with deformation-induced damage［J］.International Journal of Damage Mechanics,2006,15（2）：121- 132.

［10］MULLINS L.Softening of rubber by deformation［J］.Rubber Chemistry and Technology,1969,42（1）：339- 362.

［11］王宜舉,修乃華.非線性最優化理論與方法［M］.北京：科學出版社,2012.

［12］吳善躍,朱石堅.鋼絲繩減振器非對稱遲滯模型的參數分解識別研究［J］.振動工程學報,2007,20（06）：613- 616.

WU Shan-yue,ZHU Shi-jian.Separated identification of asymmetric hysteretic model parameters for a wirecable vibration isolator［J］.Journal of Vibration Engineering,2007,20（06）：613- 616.

Modeling and identification of asymmetric hysteresis for Kevlar tether

YI Lin,WANG Ban,GUO Ji-feng
（College of Electrical Engineering,Zhejiang University,Hangzhou 310027,China）

A generalized Bouc-Wen model with flexible shape control was proposed to describe the asymmetric hysteresis effect of Kevlar tether under quasi-static loading.An appropriate method to identify the model parameters step by step was developed by combining the mathematical model with the features of the experimental data.Theoretical curve under cyclic load was calculated using numerical method.Results showed that the experimental data fitted well with the theoretical curve,which verified the correctness of the model and the validity of the identification method.

hysteresis effect；Bouc-Wen model；Kevlar tether；identification of parameter

10.3785／j.issn.1008-973X.2015.07.024

TB 332；TB 125

A

1008- 973X（2015）07- 1376- 06

2014- 06- 27. 浙江大學學報(工學版)網址：www.journals.zju.edu.cn／eng

國家“863”高技術研究發展計劃資助項目（2013AA7044026）；國家自然科學基金資助項目（51475411）.

易琳（1986－）,男,博士生,從事空間繩系機構及控制技術的研究.E-mail：yilinyiluo＠163.com

郭吉豐,男,教授.E-mail：gjf＠zju.edu.cn