李東武 徐超
(西北工業(yè)大學(xué)航天學(xué)院, 西安 710072)
?
基于時頻域交替法的遲滯非線性振動系統(tǒng)的穩(wěn)態(tài)響應(yīng)分析?
李東武徐超?
(西北工業(yè)大學(xué)航天學(xué)院, 西安710072)
連接界面的黏滑、摩擦行為不僅是引起結(jié)構(gòu)剛度和阻尼非線性的主要原因,而且是結(jié)構(gòu)無源阻尼的主要來源.Iwan模型能夠較好地復(fù)現(xiàn)連接界面的黏滑、摩擦行為.本文采用時頻域交替法(AlternatingFrequency/TimeDomainMethod,AFT)研究含Iwan非線性模型的單自由度振子系統(tǒng)的穩(wěn)態(tài)響應(yīng).時頻域交替法具有頻域法求解線性系統(tǒng)響應(yīng)的高效性和時域法判斷非線性力的便捷性特點,采用離散傅里葉變換和傅里葉逆變換,在頻域和時域內(nèi)分別求解系統(tǒng)響應(yīng)和對應(yīng)的非線性恢復(fù)力,再反復(fù)迭代計算系統(tǒng)的穩(wěn)態(tài)響應(yīng).將時頻域交替法計算結(jié)果和中心差分法計算的結(jié)果進行對比,并研究激勵幅值對系統(tǒng)非線性特征的影響.結(jié)果表明,時頻域交替法計算的結(jié)果與中心差分計算的結(jié)果具有較好的一致性,且求解效率較高,計算耗時減少50%;隨著激勵幅值的增加,系統(tǒng)的能量耗散增加,剛度降低,固有頻率降低.
連接,遲滯非線性,Iwan模型,時頻域交替法,穩(wěn)態(tài)響應(yīng)
引言
航空、航天和機械等工程裝備系統(tǒng)中的部件大多是通過螺栓、鉚釘?shù)冗B接裝置裝配而成.當連接結(jié)構(gòu)受到外部周期載荷激勵時,連接結(jié)合面間會發(fā)生摩擦、滑移和間隙等非線性行為,這些非線性行為不僅是引起結(jié)構(gòu)振動能量耗散和阻尼衰減的主要原因,而且還造成了連接結(jié)構(gòu)剛度和阻尼的非線性.連接面間的摩擦行為是一個非常復(fù)雜的物理現(xiàn)象,它是系統(tǒng)無源阻尼的主要來源[1].已有研究表明,在切向載荷作用下,連接界面的運動將經(jīng)歷從局部微觀滑移到宏觀相對運動的過程,切向恢復(fù)力與相對位移之間為非線性遲滯關(guān)系[2-3].為了描述界面的這種多尺度粘滑摩擦過程,國內(nèi)外學(xué)者已提出一些唯象的遲滯非線性模型,例如,美國sandia國家實驗室提出的四參數(shù)Iwan模型,德國斯圖加特大學(xué)提出的Valanis模型以及微分形式的Bouc-Wen模型等[4-5].其中,Iwan模型形式簡單,模型參數(shù)具有一定的物理意義,被認為是目前研究連接界面遲滯非線性粘滑摩擦行為的一種較好的模型.
Iwan遲滯非線性模型是一種不光滑模型,因此系統(tǒng)的動力學(xué)求解較為復(fù)雜.目前常用的求解非線性系統(tǒng)穩(wěn)態(tài)響應(yīng)的方法有兩種:第一種是時域法,比如中心差分法、Runge-Kutta法等直接積分方法.直接積分法用有限元法離散空間變量,用有限差分法離散時間變量,它適合于任何情況.這類方法的優(yōu)點是只要時間步長選擇合適,計算結(jié)果會有很高的精度,而且可以運用于任意復(fù)雜的非線性系統(tǒng);缺點是在速度轉(zhuǎn)換點附近需要采用較小的步長以保證計算結(jié)果的精度,導(dǎo)致計算過程所需的時間較長.第二種是頻域法,比如諧波平衡法[6].該方法先假設(shè)出系統(tǒng)響應(yīng)的諧波解,再代入振動系統(tǒng)方程進行傅里葉級數(shù)展開求解未知參數(shù)來得到系統(tǒng)的穩(wěn)態(tài)響應(yīng).這類方法的特點是:低階諧波下計算效率較高但精度較差,隨著諧波階次的增加,計算精度變高但計算量迅速增大,計算效率降低.
時頻交替法(AlternatingFrequency/TimeDomainMethod)是一種求解復(fù)雜非線性系統(tǒng)動力響應(yīng)的新方法,其求解過程融合了頻域中求解系統(tǒng)振動微分方程的高效性(使微分方程運算變?yōu)閿?shù)值微分方程)和時域內(nèi)判斷非線性力的便捷性,并通過快速離散傅里葉正逆變換,將時域函數(shù)和頻域函數(shù)進行轉(zhuǎn)換,經(jīng)過反復(fù)的迭代最終得到振動系統(tǒng)的穩(wěn)態(tài)響應(yīng).本文采用時頻域交替法來求解含Iwan非線性模型單自由度系統(tǒng)的穩(wěn)態(tài)響應(yīng),首先給出了Iwan串-并聯(lián)模型的非線性力表達式,推導(dǎo)了時頻域交替法求解的具體公式和步驟,并對時頻域交替法求解計算的收斂性和效率進行了分析研究.
圖1所示為Iwan串-并聯(lián)模型,它是由N個Jenkins單元并聯(lián)而成,Jenkins單元則是由剛度為kn/N的線性彈簧和屈服力為fi?/N的阻尼滑塊串聯(lián)而成.在振動過程中,如果只有部分Jenkins單元發(fā)生滑移,則稱模型發(fā)生了微滑移,微滑移現(xiàn)象普遍存在,它是結(jié)構(gòu)無源阻尼的重要來源;如果所有Jenkins單元都發(fā)生滑移,則稱模型發(fā)生了宏滑移.此時模型恢復(fù)力f 達到最大值fy(稱為臨界宏觀滑移力).圖2為Iwan模型微觀滑動時的力-位移關(guān)系曲線(遲滯回線),包括初始加載段oa(又稱骨干曲線)、卸載段abc和加載段cda.
圖1 Iwan模型Fig.1 Iwanmodel圖2 遲滯曲線Fig.2 Hysteresisloop
本文選用連續(xù)Iwan模型作為研究對象,“連續(xù)”是指模型中Jenkins單元的數(shù)目N趨于無窮,采用概率密度函數(shù)φ(f?)來定義阻尼滑塊的屈服力fi?.在螺栓連接模型中,結(jié)合面接觸壓強從螺栓連接處沿連接面展向呈線性單調(diào)遞減,而且模型屈服力和接觸壓強具有相同的分布形式[7-8],因此,屈服力fi?(i=1,2,…)在[0,fN]區(qū)間內(nèi)服從均勻分布.假設(shè)Iwan模型發(fā)生宏觀滑移,則:
(1)
解得
fN=2fy
(2)
因此,根據(jù)均勻分布的定義可以得到屈服力的概率密度函數(shù)φ(f?)為
(3)
當Iwan模型受到外激勵作用時,其恢復(fù)力為
(4)
將式(3)代入式(4)得到恢復(fù)力函數(shù):
(5)
(6)
(7)
上式中左右箭頭分別表示卸載和加載,A表示模型滑動的最大位移量.
從模型的力-位移關(guān)系曲線圖2可以看出,Iwan模型的能量耗散為曲線所包絡(luò)的面積, 下面根據(jù)系統(tǒng)的穩(wěn)態(tài)響應(yīng),來推導(dǎo)單位周期內(nèi)Iwan模型的能量耗散為:
(8)
將式(6)、(7)代入式(8)得微觀黏滑下的能量耗散為:
(9)
同理可得宏觀滑移下的能量耗散為:
(10)
在結(jié)構(gòu)動力學(xué)中,頻域結(jié)果比時域結(jié)果包含更多的結(jié)構(gòu)響應(yīng)信息,在頻域中求解線性系統(tǒng)可以有效避免卷積計算[10-11],因此頻率法經(jīng)常被用來求解動力學(xué)方程.通常情況下,非線性力是系統(tǒng)狀態(tài)軌跡函數(shù),它不是時間的顯式函數(shù),這就使得非線性力很難在頻域內(nèi)分析,但是它卻容易在時域內(nèi)計算獲得.時頻域交替法融合了時域法和頻域法的優(yōu)點,利用了頻域法求解線性系統(tǒng)響應(yīng)的高效性和時域法計算非線性力的便捷性,通過快速傅里葉變換(FFT)將非線性力的時域函數(shù)轉(zhuǎn)換至頻域,求解得到系統(tǒng)的頻域穩(wěn)態(tài)響應(yīng),再將頻域響應(yīng)通過傅里葉逆變換轉(zhuǎn)換至?xí)r域,進而得到系統(tǒng)的非線性力時域函數(shù),通過反復(fù)迭代最終得到系統(tǒng)的穩(wěn)態(tài)解.在該方法計算時,需要將變量反復(fù)在時域和頻域間相互轉(zhuǎn)換.
圖3所示為單自由度振子系統(tǒng),其中m、k、c分別為系統(tǒng)振子質(zhì)量、剛度、阻尼,f(t)為外激勵,x(t)為振子的位移響應(yīng),Iwan代指圖1所示的Iwan非線性模型.
圖3 含Iwan非線性模型的單自由度振子系統(tǒng)Fig. 3 Single DOF oscillator system with Iwan nonlinear model
圖3所示單自由度系統(tǒng)的振動微分方程為:
(11)
(12)
將其代入到式(11)得到頻域內(nèi)的振動方程:
(-mω2+jωc+k)·x(ω)=
f(ω)-fn(ω,x(ω))
(13)
其中ω是指由快速傅里葉變換引入的離散頻率:
ω={ωk}
(14)
(15)
其中ΔT為采樣周期,N為快速傅里葉變換中時域函數(shù)的樣本數(shù)目,由于快速傅里葉變換后的頻域值關(guān)于ωN/2對稱,因此k的取值長度為N/2+1.
(-mω2+jωc+k)·xi(ω)=
f(ω)-fn(ω,xi-1(ω))
(16)
在這個過程中,初始位移響應(yīng)和速度響應(yīng)的獲取方法是:忽略非線性力的影響,將系統(tǒng)看作簡單的線性系統(tǒng),通過時頻域交替法求得系統(tǒng)的穩(wěn)態(tài)響應(yīng).即
(-mω2+jωc+k)·x(ω)=f(ω)
(17)
算法計算過程如圖4所示。
圖4 時頻域交替法計算的流程圖Fig.4 Flow diagram of the calculation through AFT method
算法計算過程中的關(guān)鍵控制參數(shù)有采樣頻率、決定迭代收斂與否的控制參數(shù),考慮到算法計算精度和計算效率的要求,本文選擇采樣頻率為128Hz,也可選為256Hz;迭代收斂的控制參數(shù)由第4節(jié)詳細介紹.
考慮兩種激勵幅值下的系統(tǒng)響應(yīng),第一種情況:f(t)=1.5sin(3t),時頻域交替法求解的結(jié)果如圖5、圖6所示 :
圖5 模型微滑時的穩(wěn)態(tài)響應(yīng)Fig. 5 Steady-state response in microslip
圖6 模型微滑時的遲滯曲線Fig. 6 Hysteresis loop in microslip
從圖6可以看出,在此種情形下,由于激勵幅值較小,使得非線性力的最大值沒有達到宏觀滑移力,Iwan非線性模型發(fā)生微滑移現(xiàn)象.表1列出了激勵幅值為1.5N時時頻域交替法和中心差分法所得結(jié)果的比較.
表1 模型微滑時不同方法所得結(jié)果的比較
第二種情況:f(t)=5sin(3t),結(jié)果下圖所示.
圖7 模型宏滑時的穩(wěn)態(tài)響應(yīng)Fig.7 Steady-state response in macroslip
圖8 模型宏滑時的遲滯曲線Fig. 8 Hysteresis loop in macroslip
當激勵幅值增大后,系統(tǒng)穩(wěn)態(tài)響應(yīng)相應(yīng)增大,從而使Iwan模型中更多的Jenkins單元發(fā)生滑移;從圖8可以看出,此種情形下,非線性力的最大值達到宏觀滑移力,Iwan模型發(fā)生宏觀滑移.表2相較于表1結(jié)果誤差較為明顯,但此時計算精度仍很高.誤差原因可能是非線性力的增大使得其對結(jié)果的影響程度變大.
表2 模型宏滑時不同方法所得結(jié)果的比較
運用中心差分法對時頻域交替法的計算結(jié)果進行驗證,從上面兩組計算結(jié)果可以看出,時頻域交替法求得的單自由度振子系統(tǒng)穩(wěn)態(tài)響應(yīng)能與中心差分法的結(jié)果很好地吻合,說明了時頻域交替算法適合求解此類非線性問題,而且計算的精度比較高.
下面通過求解系統(tǒng)的頻響曲線和能量耗散,進一步研究Iwan非線性模型的特性.
從圖9中可以看出,隨著激勵幅值的增大,頻響曲線的峰值逐漸向左移動,而且頻響曲線出現(xiàn)明顯的軟化現(xiàn)象,說明系統(tǒng)的固有頻率降低了.導(dǎo)致這種現(xiàn)象的本因是Iwan模型的非線性特性,具體表現(xiàn)為激勵幅值的增加使得Iwan模型中發(fā)生滑動的Jenkins單元數(shù)增加,從而導(dǎo)致振子系統(tǒng)的剛度降低,固有頻率降低.在某些參數(shù)下,該模型的頻響曲線也會發(fā)生非線性跳躍和多值現(xiàn)象.
圖9 不同激勵幅值下的頻響曲線Fig.9 Frequency response diagrams with different excitation amplitudes
隨著激勵幅值的增大到一定程度,頻響曲線中的部分頻段內(nèi)的響應(yīng)就逐漸表現(xiàn)出線性系統(tǒng)的性質(zhì),即這一段內(nèi)的不同激勵下的頻響曲線重合了.
從圖10中可以看出,Iwan非線性模型的能量耗散隨著系統(tǒng)穩(wěn)態(tài)響應(yīng)的逐漸增大而增大.微觀黏滑狀態(tài)下,非線性模型的能量耗散較小,且變化不大;宏觀滑移狀態(tài)下,非線性模型的能量耗散急劇增加.如果在雙對數(shù)坐標系下繪制能量耗散,則曲線的斜率為3,這與一些連接結(jié)構(gòu)的實驗結(jié)果是一致的.
圖10 模型能量耗散隨響應(yīng)幅值的變化曲線Fig.10 Energy dissipation-amplitude relationships
在時頻域交替法計算過程中,如何控制程序的收斂是保證算法計算效率的關(guān)鍵.本文選擇的收斂性判別要求有兩個,一是計算精度,包括兩種算法穩(wěn)態(tài)響應(yīng)的相位差和穩(wěn)態(tài)響應(yīng)的幅值比較;另一個是計算效率,指計算所需時間的比較.
選用的收斂準則為:取前后兩次迭代所得的穩(wěn)態(tài)響應(yīng)的一半長度(為了保證響應(yīng)為穩(wěn)態(tài)響應(yīng))進行相對偏差的分析,記前一次迭代得到的位移響應(yīng)為數(shù)組A,后一次迭代得到的位移響應(yīng)為數(shù)組B.利用相對偏差的大小來判斷收斂性.相對偏差的定義為
(8)
選用激勵f(t)=20sin(0.25t)來進行計算,取時間長度6000s,對不同的相對偏差值所對應(yīng)的結(jié)果進行比較.
表3 不同er的相對誤差比較
根據(jù)表3中數(shù)據(jù)的對比可以看出,相對偏差越小,時頻域交替法的計算精度越高、但是計算耗費越大;本文選擇相對偏差er=0.1作為算法迭代收斂的準則,這樣既可以滿足很高的計算精度,也能相對中心差分法而言提高計算效率.當然,如果需要保證更高的精度,可以選擇更小的收斂準則,但是計算耗費則會很大.
此外,響應(yīng)幅值和相位存在的微小誤差可以通過調(diào)節(jié)采樣頻率進行降低.
1)本文推導(dǎo)了時頻域交替法求解含Iwan非線性模型單自由度系統(tǒng)穩(wěn)態(tài)響應(yīng)的具體公式和步驟.通過文中算例的驗證及時頻域交替法收斂性和計算精度的分析,顯示出時頻域交替法計算的結(jié)果與中心差分法的結(jié)果具有很好的一致性,說明了時頻域交替法在求解非線性振動問題上具有明顯的優(yōu)勢,不僅能滿足較高的計算精度要求,而且求解效率較高,計算耗時減少約50%.
2)采用時頻域交替法求解計算,得到了單自由度振子系統(tǒng)受迫振動的幅頻響應(yīng)曲線,圖中表現(xiàn)出明顯的軟化特征,隨著激勵幅值的增加,系統(tǒng)的能量耗散增加,發(fā)生滑動的Jenkins單元數(shù)隨之增加,系統(tǒng)剛度降低,固有頻率降低.
1GaulL,LenzJ.Nonlineardynamicsofstructureassembledbyboltedjoints. Acta Mechanics Sinica,1997,125:169~181
2BergerEJ.Frictionmodelingfordynamicsystemsimulation. Journal of Applied Mechanics Reviews, 2002,55(6):535~577
3IwanWD.Adistributed-elementmodelforhysteresisanditssteay-statedynamicresponse. Journal of Applied Mechanics,1966,33(4):893~900
4ValanisKC.Atheoryofvisco-plasticitywithoutayieldsurface. Archive Mechanics, 1971,23(4):517~551
5WenYK.Methodofrandomvibrationofhysteresissystem. Archiwum Mechanics Division,1976,102(2):219~263
6SanliturkKY,ImregunM,EwinsDJ.Harmonicbalancevibrationanalysisofturbinebladeswithfrictiondampers. Journal of Vibration and Acoustics,1997,119(1):96~103
7HeK,ZhuWD.FiniteelementmodelingofstructureswithL-shapedbeamsandboltedjoints. Journal of Vibration and Acoustics-Transactions of the ASME, 2011,133(1):1~138GroperM.Microslipandmacroslipinboltedjoints. Experimental Mechanics,1985,25(2):171~174
9SegalmanDJ.ModelingjointfrictioninstructuralDyna-mics. Structural Control & Health Monitoring,2006,13(1):430~453
10CameronTM,GriffinJH.Analternatingfrequency/timedomainmethodforcalculatingthesteady-stateresponseofnonlinearsystems. Journal of Applied Mechanics,1989,56(3):149~154
11朱孟華. 求解非線性系統(tǒng)振動響應(yīng)的Fourier正逆變換方法. 內(nèi)燃機工程,1995,16(3):1~8 (ZhuMH.Amethodforsolvingthenon-linearvibrationsystemsbyfourierdirectinversetransformtechnique. Chinese Internal Combustion Engine Engineering, 1995,16(3):1~8 (inChinese))
*TheprojectsupportedbytheNationalNaturalScienceFoundationofChina(11372246).
?CorrespondingauthorE-mail:chao_xu@nwpu.edu.cn
02June2014,revised29June2014.
ALTERNATINGTIME/FREQUENCYDOMAINMETHODFORCALCULATINGTHESTEADY-STATERESPONSEOFHYSTERESISNONLINEARVIBRATIONSYSTEMS?
LiDongwuXuChao?
(School of Astronautics, Northwestern Polytechnical University, ShanXi, Xian710072, China)
Stick-slipandfrictionofthejointinterfacepayasignificantrolenotonlyonstructurenonlinearstiffnessanddampingbutalsoonstructurepassivedamping.Iwanmodelprovidesapreferablerecurrenceofthestick-slipandfrictioncharacteristicsofthejointinterface.Thealternatingfrequency/timedomainmethod(AFT)isusedtostudythesteady-stateresponseofsingledegreeoffreedomoscillatorsystemcontainingIwannonlinearmodelinthispaper.Thismethodtakesbothadvantagesoftheeffectivenessofcalculatingresponseforlinearsystembyfrequencydomainmethodandtheeaseofevaluatingnonlinearforcebytimedomainmethod.ThediscreteFouriertransformalternatingfromthetimedomaintothefrequencydomainisappliedthroughiteratingtoobtainthesteady-stateresponseofthesystem.TheresultofAFTiscomparedwiththatfromthecentraldifferencemethod,andtheinfluenceofexcitationamplitudeonnonlinearcharacteristicisalsoexaminedinthispaper.Theresultsshowthatthesteady-stateresponsebyAFTmethodagreeswellwiththatofthecentraldifferencemethod.Moreover,theAFTmethodperformsbettercomputationalefficiency,asthetime-consumingisreducedby50percent.Inawhole,itisfoundthattheapplicationofhigherexcitationamplituderesultsintheincreaseofthesystemenergydissipationandthereducesstiffnessandnaturalfrequency.
joint,hysteresisnonlinear,Iwanmodel,alternatingfrequency/timedomainmethod,steady-stateresponse
E-mail:chao_xu@nwpu.edu.cn
10.6052/1672-6553-2015-043
2015-06-02收到第1稿,2015-06-29收到修改稿.
*國家自然科學(xué)基金資助項目(11372246)