閻 俊，徐 超

(西北工業(yè)大學(xué) 航天學(xué)院，西安 710072)

工程結(jié)構(gòu)多由螺栓、鉚釘?shù)冗B接裝配。振動(dòng)環(huán)境下連接結(jié)合面在法向可能出現(xiàn)間隙、分離及碰撞，在切向可能出現(xiàn)摩擦、滑移及滑動(dòng)等非線性現(xiàn)象，導(dǎo)致連接結(jié)構(gòu)發(fā)生復(fù)雜的動(dòng)力學(xué)行為。對(duì)含結(jié)合面摩擦系統(tǒng)的研究一直為非線性動(dòng)力學(xué)領(lǐng)域研究重點(diǎn)。何阿平等[1]對(duì)摩擦學(xué)行為的混沌特性進(jìn)行定量研究，并實(shí)驗(yàn)證明摩擦學(xué)系統(tǒng)混沌特性。吳文健等[2]分析由干摩擦引起的粘滑振動(dòng)及系統(tǒng)穩(wěn)定性，采用數(shù)值方法獲得系統(tǒng)經(jīng)周期運(yùn)動(dòng)失去穩(wěn)定性通用混沌路徑。韓彥等[3]對(duì)干摩擦振動(dòng)系統(tǒng)響應(yīng)計(jì)算方法進(jìn)行述評(píng)。陳建縣[4]分別從摩擦特性、摩擦力模型、摩擦系統(tǒng)自激振動(dòng)、強(qiáng)迫振動(dòng)及摩擦振動(dòng)控制等系統(tǒng)綜述機(jī)械系統(tǒng)滑動(dòng)干摩擦動(dòng)力學(xué)近期研究進(jìn)展，指出含摩擦系統(tǒng)動(dòng)力學(xué)研究對(duì)工程設(shè)計(jì)、優(yōu)化及控制的重要意義。

結(jié)合面干摩擦描述形式對(duì)建模系統(tǒng)非線性動(dòng)力學(xué)特性有重要影響[5]。文獻(xiàn)[6]對(duì)機(jī)械系統(tǒng)中摩擦模型進(jìn)行評(píng)述，系統(tǒng)介紹多種常用摩擦模型，并對(duì)每種模型的構(gòu)成、特點(diǎn)及適用范圍進(jìn)行較詳細(xì)論述，認(rèn)為動(dòng)態(tài)摩擦模型能更好描述摩擦非線性行為。文獻(xiàn)[7]對(duì)采用不同動(dòng)態(tài)摩擦模型的系統(tǒng)非線性動(dòng)力學(xué)行為進(jìn)行研究，指出不同摩擦力模型的優(yōu)點(diǎn)與不足。對(duì)一般機(jī)械結(jié)合面而言，接觸壓力分布不均勻，在切向振蕩載荷作用下，結(jié)合面將反復(fù)經(jīng)歷從微觀尺度滑移到宏觀尺度滑動(dòng)的多尺度動(dòng)態(tài)摩擦過(guò)程。機(jī)械摩擦動(dòng)力學(xué)研究中廣泛采用的庫(kù)倫摩擦模型雖形式簡(jiǎn)單，卻不能準(zhǔn)確描述多尺度物理過(guò)程[4,6]。文獻(xiàn)[8]從結(jié)合面摩擦物理機(jī)理入手，在遲滯非線性模型[9]基礎(chǔ)上導(dǎo)出描述連接面多尺度粘滑摩擦過(guò)程的參數(shù)化模型形式。文獻(xiàn)[10]通過(guò)引入并聯(lián)線性剛度項(xiàng)及均勻分布臨界滑移力分布密度函數(shù)形式，給出改進(jìn)的描述結(jié)合面多尺度粘滑過(guò)程的力學(xué)模型。

本文以典型含連接結(jié)合面振動(dòng)系統(tǒng)為研究對(duì)象，采用滑動(dòng)位移的遲滯模型[10]描述結(jié)合面多尺度粘滑摩擦過(guò)程，利用中心差分法，數(shù)值研究系統(tǒng)受諧波強(qiáng)迫激勵(lì)時(shí)在主共振、次共振條件下的混沌特性。

1 含結(jié)合面多尺度粘滑摩擦系統(tǒng)模型

機(jī)械結(jié)構(gòu)中典型的螺栓連接見(jiàn)圖1。兩被連接件通過(guò)施加預(yù)緊載荷螺栓裝配在一起。通常預(yù)緊力在結(jié)合面上產(chǎn)生的接觸壓力p分布并非均勻，靠近螺栓孔中心區(qū)域的接觸壓力高，離中心越遠(yuǎn)接觸壓力越低。

圖1 切向力作用螺栓連接結(jié)構(gòu)

受切向載荷作用時(shí)，接觸壓力較低區(qū)域先發(fā)生滑動(dòng)，而接觸壓力較高區(qū)域仍保持粘滯，從而接觸界面處于既存在局部滑移亦存在局部粘滯狀態(tài)；若切向載荷足夠大，則整個(gè)結(jié)合面處處均會(huì)發(fā)生滑動(dòng)，即形成動(dòng)態(tài)相對(duì)運(yùn)動(dòng)。局部滑移可能發(fā)生于10-6～10-4m尺度，而宏觀滑動(dòng)則可能發(fā)生在10-3～10-2m尺度；因此，結(jié)合面動(dòng)態(tài)摩擦過(guò)程具有典型的多尺度特征。

考慮圖1諧波激勵(lì)下典型的螺栓連接系統(tǒng)，因?qū)ΨQ性且只考慮沿結(jié)合面的切向運(yùn)動(dòng)，可簡(jiǎn)化為圖2(a)的運(yùn)動(dòng)模型。模型中下連件固定，上連件考慮為一維運(yùn)動(dòng)剛體。u，f，fesin(ωt)分別為結(jié)合面相對(duì)位移、非線性遲滯恢復(fù)力及諧波外激勵(lì)力。連接系統(tǒng)剛度包括并聯(lián)兩部分，即由結(jié)合面粘滑摩擦行為的貢獻(xiàn)及因連接件柔性的貢獻(xiàn)。前者隨結(jié)合面相對(duì)位移增大而減小呈非線性變化，后者視作線性剛度。令k表示系統(tǒng)初始總剛度，與結(jié)合面粘滑摩擦行為無(wú)關(guān)的線性剛度為ka=k，無(wú)量綱系數(shù)0<α<1；則描述結(jié)合面粘滑摩擦部分的子系統(tǒng)初始剛度為k-ka。

圖2 Iwan摩擦振子及Iwan模型

力學(xué)模型實(shí)際為并聯(lián)的線性剛度項(xiàng)Iwan模型[10]，見(jiàn)圖2(b)，其采用n個(gè)彈簧滑塊單元并聯(lián)組成的子系統(tǒng)(經(jīng)典Iwan模型)描述結(jié)合面多尺度粘滑摩擦行為。子系統(tǒng)中每個(gè)線性彈簧剛度ki(i=1,2,…n)均為(k-ka)/n，但每個(gè)滑塊的臨界滑移力fi(i=1,2,…n)不相同。系統(tǒng)切向受載時(shí)，臨界滑移力小的滑塊先發(fā)生滑動(dòng)，隨相對(duì)位移增大，發(fā)生滑移滑塊越多，直至全部滑移，即結(jié)合面出現(xiàn)宏觀相對(duì)運(yùn)動(dòng)；因此，該模型可較好復(fù)現(xiàn)結(jié)合面上復(fù)雜的多尺度粘滑摩擦過(guò)程。

考慮經(jīng)典Iwan模型連續(xù)形式，取n→∞，則fi需定義成分布密度函數(shù)Φ(fi)形式(表達(dá)式非唯一)。定義Φ(fi)為均勻分布形式[11]，即

(1)

式中：fy為結(jié)合面恰發(fā)生宏觀滑移時(shí)對(duì)應(yīng)的臨界宏觀滑移力。

經(jīng)典Iwan模型恢復(fù)力骨干(backbone)曲線fo可表示為未發(fā)生滑移的彈簧滑塊單元與發(fā)生滑移的彈簧滑塊單元上恢復(fù)力之和,即

(2)

將式(1)代入式(2)積分并整理得：

fo[k(1-α),u]=

(3)

式(3)分段函數(shù)包括兩部分，分別表示結(jié)合面處于局部滑移與宏觀滑動(dòng)狀態(tài)時(shí)恢復(fù)力-變形關(guān)系。振蕩載荷下結(jié)合面會(huì)出現(xiàn)反復(fù)的卸載-加載運(yùn)動(dòng)。據(jù)梅辛準(zhǔn)則[8],得遲滯過(guò)程中Iwan模型恢復(fù)力函數(shù)為

(4)

式中：fl，fu分別為卸載、加載中恢復(fù)力；A為系統(tǒng)一個(gè)振動(dòng)周期內(nèi)最大振幅。

引入無(wú)量綱參數(shù)g=f/fy,q=ku/fy,Q=kA/fy對(duì)式(2)、(4)正則化得：

(5)

(6)

據(jù)牛頓第二運(yùn)動(dòng)定律，圖2(a)系統(tǒng)運(yùn)動(dòng)方程可表示為

(7)

對(duì)式(7)正則化，兩邊同乘k/(mfy)得：

(8)

(9)

2 計(jì)算方法

式(9)振動(dòng)系統(tǒng)遲滯恢復(fù)力具有變剛度、變阻尼非光滑、強(qiáng)非線性特點(diǎn)，公式解析解較難獲得。本文用中心差分法求解運(yùn)動(dòng)方程，研究多尺度粘滑摩擦系統(tǒng)隨外激勵(lì)幅值逐漸增大表現(xiàn)出的分岔及混沌特征。

由于恢復(fù)力分段及遲滯特性，數(shù)值仿真中準(zhǔn)確判斷轉(zhuǎn)折點(diǎn)至關(guān)重要。

2.1 局部滑移與宏觀滑移轉(zhuǎn)折點(diǎn)

由式(5)知，該轉(zhuǎn)折點(diǎn)對(duì)應(yīng)的結(jié)合面相對(duì)位移為

(10)

當(dāng)結(jié)合面相對(duì)位移超過(guò)qmax時(shí)，則判斷系統(tǒng)進(jìn)入宏觀滑移狀態(tài)，恢復(fù)力用式(5)第2式計(jì)算。

2.2 遲滯回線卸-加載轉(zhuǎn)折點(diǎn)

遲滯回線中卸-加載轉(zhuǎn)折點(diǎn)即系統(tǒng)相對(duì)速度為零時(shí)刻，判斷式為

(11)

即當(dāng)i與i+1時(shí)刻速度符號(hào)相反、兩者乘積滿足小于等于零，且要求i+1時(shí)刻速度不等于零時(shí)，則i時(shí)刻為卸-加載轉(zhuǎn)折點(diǎn)。由于數(shù)值仿真中采用固定步長(zhǎng)，式(11)判斷的轉(zhuǎn)折點(diǎn)僅為近似轉(zhuǎn)折點(diǎn)，但系統(tǒng)分析步長(zhǎng)足夠小時(shí)，近似誤差亦小。

(12)

當(dāng)ε=0時(shí)，式(12)線性系統(tǒng)稱為式(7)的派生系統(tǒng)，其中ωo為派生系統(tǒng)固有頻率，ω為外激勵(lì)頻率。ωo接近外激勵(lì)頻率ω時(shí)(ωo≈ω)，系統(tǒng)存在主共振現(xiàn)象；ωo接近外激勵(lì)頻率ω的1/2時(shí)(2ωo≈ω)，系統(tǒng)發(fā)生強(qiáng)烈1/2次亞諧波共振現(xiàn)象。引入比值變量Ω=ω/ωo，取Ω=1，Ω=2分別研究系統(tǒng)主共振及1/2亞諧波共振下混沌形態(tài)。每類(lèi)情況下分別選α=0(不考慮線性剛度項(xiàng))與α=0.5兩種，對(duì)比分析發(fā)生宏觀滑移后系統(tǒng)剩余剛度對(duì)混沌形態(tài)影響。

令時(shí)間τ為第三維，將式(9)平面非自治非線性系統(tǒng)擴(kuò)維成三維，得系統(tǒng)龐加萊截面為

(13)

建立τ=0時(shí)初值點(diǎn)及其按解的路徑經(jīng)時(shí)間2π/Ω后所達(dá)點(diǎn)間之關(guān)系。即主共振與1/2亞諧共振下系統(tǒng)平衡點(diǎn)附近的龐加萊映射。進(jìn)一步將系統(tǒng)離散成龐加萊點(diǎn)映射系統(tǒng)為

gi+1=y(gi)

(14)

(15)

3 結(jié)果與討論

3.1 1∶1共振

α=0時(shí)，圖3(a)～(d)給出典型參數(shù)值下系統(tǒng)龐加萊截面。全局分岔見(jiàn)圖7(a)。由圖3可知，ge在零點(diǎn)附近激勵(lì)幅值量級(jí)較小不足以使滑塊單元產(chǎn)生滑動(dòng)，系統(tǒng)非線性不明顯。在微小擾動(dòng)下僅存在s+，s-兩段實(shí)時(shí)變化的互異雙曲型平衡曲線(圖3(a) 龐加萊映射)反對(duì)稱落在q=0軸兩邊；ge繼續(xù)增大，龐加萊截面上兩曲線段在ge=2處螺旋式趨于各自與之同側(cè)的互異雙曲平衡點(diǎn)s+，s-(圖3(b)龐加萊映射)，此處兩平衡點(diǎn)穩(wěn)定流形Ws(s+)，Ws(s-)與不穩(wěn)定流形Wu(s+)，Wu(s-) 滿足：Ws(s+)∩Wu(s-)≠0，或Ws(s-)∩Wu(s+)≠0，即重合后形成不變流形，由此斷定ge=2為異宿點(diǎn)；異宿軌道上點(diǎn)當(dāng)τ→+∞及τ→-∞時(shí)分別趨近于各自不同的s+，s-。上述流形穩(wěn)定性判別可據(jù)Peixoto[12]定理；由于該系統(tǒng)不同雙曲鞍點(diǎn)s+，s-的穩(wěn)定流形Ws(s+)，不穩(wěn)定流形Wu(s-)橫截相交于ge=2點(diǎn)，則該點(diǎn)亦稱為橫截異宿點(diǎn)，可通過(guò)Melnikov[13]方法判定。ge逐漸增大到2.07，龐加萊截面上平衡點(diǎn)(周期解)個(gè)數(shù)突增至6個(gè)，圖7(a)拓?fù)浞种c圖3(c)說(shuō)明經(jīng)叉式分岔打破原穩(wěn)定平衡解。ge=2.5時(shí)，系統(tǒng)由P6周期軌道過(guò)渡到混沌帶，圖3(d)～(f)分別為該點(diǎn)的龐加萊映射、功率譜、吸引子；之后系統(tǒng)重復(fù)經(jīng)歷P2→混沌過(guò)程；可見(jiàn)，系統(tǒng)伴隨剛度減小及彈塑性過(guò)程的不斷交替，最終出現(xiàn)混沌運(yùn)動(dòng)原因?yàn)闄M截異宿點(diǎn)的產(chǎn)生，此時(shí)雙曲鞍點(diǎn)穩(wěn)定流形與不穩(wěn)定流形為不同吸引子的吸引域邊界，與不穩(wěn)定流形橫截相交后變?yōu)闃O復(fù)雜從而成為分形邊界，使運(yùn)動(dòng)具有極端敏感性。

圖3 龐加萊截面、功率譜及吸引子

圖4 龐加萊截面

α= 0.5時(shí)，與α=0對(duì)比發(fā)現(xiàn)，最初小擾動(dòng)下系統(tǒng)仍存在s+，s-兩個(gè)實(shí)時(shí)變化的互異雙曲型平衡曲線段見(jiàn)圖4(a)，反對(duì)稱落在q=0軸兩邊，較α=0早于ge=1處開(kāi)始螺旋式趨于與之同側(cè)的互異雙曲平衡點(diǎn)s+，s-，見(jiàn)圖4(b)。同理可驗(yàn)證ge=1為異宿點(diǎn)，但異宿軌道上點(diǎn)τ→-∞時(shí)，趨近于s+，s-；而τ→+∞時(shí)，軌線上點(diǎn)不再趨近于s+，s-，而發(fā)散開(kāi)形成多倍周期Pn軌道，龐加萊截面上呈現(xiàn)多個(gè)離散點(diǎn)見(jiàn)圖4(c)～(d)；由此判斷ge=1非橫截異宿點(diǎn)。全局分岔見(jiàn)圖7(b)，系統(tǒng)始終未出現(xiàn)混沌軌道；因此，線性剛度項(xiàng)對(duì)系統(tǒng)動(dòng)力學(xué)行為影響至關(guān)重要。

3.2 1∶2共振

α= 0時(shí)，微小擾動(dòng)下零點(diǎn)附近存在s1，s2不變流形雙曲平衡點(diǎn)見(jiàn)圖5(a)。ge漸增至0.5，圖5(b)為兩個(gè)雙曲平衡點(diǎn)融合為s過(guò)程。此處穩(wěn)定流形Ws(s)與不穩(wěn)定流形Wu(s)滿足相交條件Ws(s) ∩Wu(s) ≠ 0直到ge≈0.8附近，螺旋式趨于完全融合的雙曲平衡點(diǎn)s，則上述集合中同宿點(diǎn)在ge≈0.8附近，此時(shí)對(duì)應(yīng)的Ws(s)與Wu(s)彼此重合。同宿軌道上點(diǎn)τ→±∞時(shí)均趨近于一點(diǎn)s。

圖5 龐加萊截面、功率譜及吸引子

圖6 龐加萊截面、功率譜及吸引子

圖7(c)為區(qū)間[0,2]上分岔圖。該區(qū)間同宿分支點(diǎn)不止一個(gè)；因此推斷[13]，若存在同宿點(diǎn)s∈Ws(s) ∩Wu(s)，則s同時(shí)在Ws(s)，Wu(s)上。由于Ws(s)，Wu(s)均為不變流形，故對(duì)整數(shù)m,gm(s)也在Ws(s)，Wu(s)上，gm(s)為同宿點(diǎn)，必有無(wú)窮多同宿點(diǎn)存在，由此證明ge∈[0, 2]區(qū)間上不止一個(gè)同宿點(diǎn)，使該區(qū)間上不變流形異常復(fù)雜。ge=2.55時(shí)，龐加萊截面上平衡點(diǎn)個(gè)數(shù)突增至3個(gè)，系統(tǒng)經(jīng)局部叉式分岔改變其拓?fù)浣Y(jié)構(gòu)，打破原穩(wěn)定平衡態(tài)，變化過(guò)程見(jiàn)圖5(c)；之后循環(huán)出現(xiàn)3次叉式分岔，分別為ge=3.9處P3→P1，ge=4.1處P1→P3，ge=4.5處P3→P1。ge相對(duì)較大時(shí)系統(tǒng)出現(xiàn)柱狀混沌域。圖5(d)～(f)分別為ge=14時(shí)的龐加萊映射、功率譜及奇異吸引子。

α= 0.5時(shí)，與α= 0 對(duì)比發(fā)現(xiàn)，初始階段小擾動(dòng)下龐加萊截面顯示準(zhǔn)周期運(yùn)動(dòng)軌道見(jiàn)圖6(a)，ge=0.3處轉(zhuǎn)換為圖6(b)的穩(wěn)定周期運(yùn)動(dòng)s，其穩(wěn)定流形Ws(s)與不穩(wěn)定流形Wu(s)彼此重合。同理ge=0.3點(diǎn)為同宿分岔點(diǎn)，s軌為同宿分岔軌道。同宿軌道上點(diǎn)τ→ ±∞時(shí)均趨近于s軌。圖7(d)說(shuō)明在區(qū)間ge∈[0, 16]上不止出現(xiàn)一個(gè)ge=0.3的同宿分岔點(diǎn)。圖6(c)ge=0.5也為同宿點(diǎn)，說(shuō)明α= 0中必存在無(wú)窮個(gè)同宿點(diǎn)使不變流形復(fù)雜程度提高，最終出現(xiàn)內(nèi)嵌式混沌軌道?？赏ㄟ^(guò)圖6(d)～(f)在ge=10處的龐加萊映射、功率譜、奇異吸引子說(shuō)明混沌現(xiàn)象。

圖7 全局分岔圖

4 結(jié) 論

以往簡(jiǎn)單的干摩擦庫(kù)倫模型振子大多驗(yàn)證通過(guò)周期倍化分岔、Hopf分岔及周期運(yùn)動(dòng)失穩(wěn)等常見(jiàn)路徑通向混沌。本文在干摩擦振子基礎(chǔ)上，通過(guò)對(duì)多尺度粘滑摩擦振子混沌分析，結(jié)論如下：

(1) 考慮結(jié)合面的多尺度粘滑干摩擦情況存在的異宿、同宿分岔系統(tǒng)運(yùn)動(dòng)具有極端敏感性，小擾動(dòng)下即可能出現(xiàn)混沌，從而驗(yàn)證同宿、異宿軌道通向混沌的路徑，并表明龐加萊映射的雙曲鞍點(diǎn)穩(wěn)定流形與不穩(wěn)定流形是否相交與混沌運(yùn)動(dòng)是否發(fā)生密切相關(guān)。

(2) 系統(tǒng)在1∶1共振下因存在橫截異宿點(diǎn)，雙曲鞍點(diǎn)穩(wěn)定流形、不穩(wěn)定流形為不同吸引子的吸引域邊界，與不穩(wěn)定流形橫截相交使系統(tǒng)拓?fù)浣Y(jié)構(gòu)異常復(fù)雜成為分形邊界，使運(yùn)動(dòng)具有極端敏感性，出現(xiàn)混沌。

(3) 系統(tǒng)1∶2共振下因存在無(wú)窮個(gè)同宿分支點(diǎn)，使不變流形異常復(fù)雜，系統(tǒng)受小擾動(dòng)時(shí)，可能出現(xiàn)混沌運(yùn)動(dòng)。

[1] 朱華, 葛世榮.摩擦學(xué)系統(tǒng)的混沌特征[J].機(jī)械工程學(xué)報(bào), 2004, 40(12):10-13.

ZHU Hua, GE Shi-rong.Chaotic characteristics of tribological systems[J].Chinese Journal of Mechanical Engineering, 2004, 40(12):10-13.

[2] 丁旺才, 張有強(qiáng), 張慶爽.含干摩擦振動(dòng)系統(tǒng)的非線性動(dòng)力學(xué)分析[J].工程力學(xué), 2008, 25(10):212-217.

DING Wang-cai, ZHANG You-qiang, ZHANG Qing-shuang.Nonlinear dynamics analysis of vibrate system with dry friction[J].EngineeringMechanics, 2008, 25(10):212-217.

[3] 白鴻柏, 黃協(xié)清.干摩擦振動(dòng)系統(tǒng)響應(yīng)計(jì)算方法研究綜述[J].力學(xué)進(jìn)展, 2001, 34(4):527-534.

BAI Hong-bai, HUANG Xie-qing.An overview on study of methods of response computation for the dry frictionally damped vibration systems[J].Advances in Mechanics, 2001, 34(4):527-534.

[4] 丁千, 翟紅梅.機(jī)械系統(tǒng)摩擦動(dòng)力學(xué)研究進(jìn)展[J].力學(xué)進(jìn)展, 2013, 43(1):112-131.

DING Qian, ZHAI Hong-mei.The advance in researches of friction dynamics in mechanics system[J].Advances in Mechanics, 2013, 43(1):112-131.

[5] 蔡力鋼，王鋒，李玲.栓接結(jié)合部動(dòng)態(tài)特性研究進(jìn)展[J].機(jī)械工程學(xué)報(bào), 2013, 49(9): 158-168.

CAI Li-gang, WANG Feng, LI Ling.Review on dynamic properties of bolted joints[J].Chinese Journal of Mechanical Engineering, 2013, 49(9): 158-168.

[6] 劉麗蘭, 劉宏昭, 吳子英,等.機(jī)械系統(tǒng)中摩擦模型的研究進(jìn)展[J].力學(xué)進(jìn)展, 2008, 38(2):201-213.

LIU Li-lan, LIU Hong-zhao, WU Zi-ying,et al.An overview of friction models in mechanical systems[J].Advances in Mechanics, 2008, 38(2):201-213.

[7]Awrejcewicz J, Olejnik P.Analysis of dynamic systems with various friction laws[J].Applied Mechanics Reviews, 2005, 58:389-410.

[8] Segalman D J.Modeling joint friction in structural dynamics[J].Structural Control & Health Monitoring, 2006, 13(1): 430-453.

[9] Iwan W D.Thedynamic response of the 1-DOF bilinear hysteretic system[J].In: Proc.3rd World Conf.Earthquake Eng., New Zealand: University of Auckland, 1965, 2: 783- 796.

[10] Song Y.Modeling, identification and simulation of dynamics of structures with joints and interfaces[D].Urbana, Illinois: University of Illinois, 2004.

[11] Iwan, W D.Adistributed-element model for hysteresis and its steady-state dynamics systems[J].ASME Journal of Applied Mechanics, 1966, 33:893-900.

[12] Peixoto M C, Peixoto M M.Structural stability in the plane with enlarged boundary conditions[J].Acad., Brazil.Sci., 1959, 31:135-160.

[13] 楊紹普, 申永軍.滯后非線性系統(tǒng)的分岔與奇異性[M].北京: 科學(xué)出版社,2003.