馬顏顏 寧麗娟

（陜西師范大學數學與信息科學學院，西安 710119）

高斯白噪聲激勵下分數階Duffing-Van der Pol系統(tǒng)的穩(wěn)態(tài)響應?

馬顏顏 寧麗娟?

（陜西師范大學數學與信息科學學院，西安 710119）

應用隨機平均法研究了高斯白噪聲激勵下含有分數階阻尼項的Duffing-Van der Pol系統(tǒng)的穩(wěn)態(tài)響應.首先應用基于廣義諧和函數的隨機平均法得到系統(tǒng)關于幅值的平均伊藤微分方程并建立相應的平穩(wěn)FPK方程，求解該平穩(wěn)FPK方程的近似理論解得到系統(tǒng)幅值的穩(wěn)態(tài)概率密度.分析幅值、位移和速度的穩(wěn)態(tài)概率密度探究分數階阻尼項以及其它參數對系統(tǒng)穩(wěn)態(tài)響應的影響.發(fā)現降低分數階的階數可以增強系統(tǒng)的響應而增大分數階的系數可以減弱系統(tǒng)響應.最后對原系統(tǒng)進行Monte Carlo數值模擬驗證近似理論解的有效性.

響應， 分數階， Duffing-Van der Pol， 高斯白噪聲， 隨機平均法， Monte Carlo模擬

引言

Einstein在1905年首次提出了隨機激勵下的強非線性系統(tǒng).它可以廣泛地應用于力學、控制、生物、航天科學、金融、社會科學以及粘彈性材料等眾多領域中［1-5］.非線性系統(tǒng)自被提出以來一直都是非常具有吸引力和發(fā)展前景的研究領域.然而隨著深入的研究，人們發(fā)現許多實際系統(tǒng)中存在記憶因子而不能用傳統(tǒng)的整數階模型進行模擬［6，7］.為解決這個問題Bagley和Torvik在文獻［8，9］中證明了分數階導數模型能夠更好地模擬粘彈性材料.另外越來越多的研究［10-14］表明含有分數階導數的非線性系統(tǒng)比傳統(tǒng)的整數階系統(tǒng)能更精確地模擬實際物理系統(tǒng)的內部性質且給出更合理的解釋.

事實上，近年來研究者們已應用多種方法研究了在不同噪聲激勵下的分數階非線性系統(tǒng)的動力學行為，取得了豐碩的成果:文獻［15］應用廣義諧波技術和能量包線隨機平均法研究了在高斯白噪聲激勵下含有分數階的Duffing振子的穩(wěn)態(tài)響應.文獻［16］討論了高斯白噪聲激勵下含有Caputo型分數階阻尼的非線性碰撞振動系統(tǒng)的隨機分岔.該文獻主要應用了Zhuravlev非光滑變換和隨機平均法.文獻［17］利用隨機平均法探究了分數階Duffing振子在諧和力和白噪聲共同激勵下的隨機跳躍和隨機分岔.應用隨機平均法文獻［18］討論了分數階單自由度非線性系統(tǒng)在高斯白噪聲激勵下的響應和穩(wěn)態(tài).在文獻［19］中應用隨機平均法研究了一系列含有Caputo型分數階導數的自激系統(tǒng)在高斯白噪聲激勵下的隨機響應.文獻［20］研究了實噪聲激勵下含有分數階阻尼項的多自由度非線性振子的首次穿越損壞.上述的研究基本上是應用標準隨機平均法或者能量包線隨機平均法.可以看到隨機平均法不僅可以簡化運動方程還可以降低方程的維數，因此經常被用來求解含有分數階導數的非線性系統(tǒng)的近似理論解.

除了隨機平均法還有許多其它有效的方法來求解分數階非線性系統(tǒng)的理論解.比如文獻［21］應用攝動法分析了含有分數階阻尼項的擬線性系統(tǒng)在泊松白噪聲激勵下的響應.文獻［22，23］應用L-P方法和多尺度方法研究含有分數階阻尼項的Duffing振子的動力學響應.頻域方法被應用于研究分數階隨機系統(tǒng)［24］.在文獻［25］中提出用特征向量展開方法和拉普拉斯變換方法研究含有分數階阻尼的隨機動力學系統(tǒng).

Duffing-Van der Pol系統(tǒng)是非常具有代表性的強非線性模型而且被廣泛應用于科學、生物、工程等很多領域.在文獻［26］中用隨機雙Duffing-Van der Pol網絡模型模擬了腦電圖信號的輸出.前面的文獻中涉及到多種噪聲，其中高斯白噪可用具體數學表達式表述，便于推導分析和運算.在這篇文章中應用隨機平均法研究高斯白噪聲激勵下分數階Duffing-Van der Pol系統(tǒng)的穩(wěn)態(tài)響應.

文章結構如下:第1節(jié)介紹分數階Duffing-Van der Pol系統(tǒng)模型和Riemann-Liouvile型分數階導數的定義.第2節(jié)應用隨機平均法求得原系統(tǒng)幅值的平均伊藤微分方程并建立相應的FPK方程，求解FPK方程得到幅值的穩(wěn)態(tài)概率密度.第3節(jié)通過穩(wěn)態(tài)概率密度分析參數變化對系統(tǒng)穩(wěn)態(tài)響應的影響，并對原系統(tǒng)進行Monte Carlo數值模擬檢驗理論結果的有效性.第4節(jié)給出文章結論.

1 系統(tǒng)模型

考慮如下高斯白噪聲激勵下含有分數階導數的Duffing-Van der Pol系統(tǒng):

這里δ1表示線性阻尼的系數；δ2為非線性阻尼的系數；ω無阻尼系統(tǒng)的自然頻率；μ為非線性的強度；g（X）＝ω2X+μX3代表回復力是強非線性函數；分數階導數DαX（t）采用Riemann-Liouville定義:

這里 α（0＜α＜1）是分數階的階數，Γ（·）是伽馬函數，χ表示分數階導數項的系數．ξ（t）為高斯白噪聲，均值和方差分別為:

其中δ（τ）是狄拉克函數．

2 隨機平均法和近似解

根據廣義諧和函數系統(tǒng)（1）在原點附近的運動可以近似為周期的，因此有如下變換［27］:

這里 U（X）是勢函數；cosΘ（t）和 sinΘ（t）表示廣義諧和函數；ν（A，Θ）和 Θ（t）分別表示系統(tǒng)的瞬時頻率和瞬時相位．

現把ν（A，Θ）展開成傅立葉級數的形式如下:

方程（10）關于 Θ 在［0，2π］上積分得到近似平均頻率，即:

把方程（12）代入方程（6）可得到 Θ（t）的近似表達式:

由于A（t）和 Γ（t）是關于時間 t的慢變過程，所以可由方程（13）得到如下近似表達式:

由方程（15）可得到 Stratonovich隨機微分方程，加上 Wong-Zakai［28］修正項并平均后可得對應的伊藤平均微分方程:

其中漂移項和擴散項分別為:

這里〈·〉Θ是關于 Θ 在［0，2π］上的平均．

考慮近似關系式（14），含有分數階導數項的〈M11〉Θ可進一步化簡為:

這里FPK方程關于A的邊界條件為當A＝0時 p＝finite以及當 A→∞時 p、?p／?A→0．依據這些邊界條件求解FPK方程（24）得到原系統(tǒng)（1）幅值的穩(wěn)態(tài)概率密度:

其中C 表示歸一化常數，m（s）和 σ2（s）的表達式由（22）式和（23）式給出．

總能H的穩(wěn)態(tài)概率密度為:

A＝U-1（H）表示 H＝U（A）的反函數；H 表示系統(tǒng)的總能．

廣義位移X和速度˙X的聯(lián)合穩(wěn)態(tài)概率密度為:

3 穩(wěn)態(tài)響應和參數分析

分析系統(tǒng)參數和噪聲強度的變化對系統(tǒng)穩(wěn)態(tài)響應的影響，并對原系統(tǒng)進行Monte Carlo數值模擬驗證理論方法的有效性［29，30］．圖 1 給出了幅值A、位移X和速度˙X（Y＝˙X）的穩(wěn)態(tài)概率密度曲線．其中實線表示理論結果，圓點表示原系統(tǒng)的Monte Carlo數值模擬結果．參數取值如下:α＝0．7，χ＝4，ω＝2，μ＝1．6， D＝0．06，δ1＝ 0．05，δ2＝-0．05．由圖 1知近似解析解與原系統(tǒng)的Monte Carlo數值模擬結果相吻合．這表明隨機平均法可以有效地求解高斯白噪聲激勵下分數階Duffing-Van der Pol系統(tǒng)的穩(wěn)態(tài)響應．

圖1（a） 幅值A的穩(wěn)態(tài)概率密度曲線Fig.1（a） Stationary probability density of amplitude A

圖1（b） 速度Y的穩(wěn)態(tài)概率密度曲線Fig.1（b） Stationary probability density of velocity Y

圖1（c） 位移X的穩(wěn)態(tài)概率密度曲線Fig.1（c） Stationary probability density of displacement X

圖2給出了對于不同分數階α的幅值A的穩(wěn)態(tài)概率密度曲線．參數取值如下:χ＝3，ω＝2，μ＝1．5，D＝0．06，δ1＝0．05，δ2＝ -0．05．其中實線表示理論結果，圓點表示原系統(tǒng)的Monte Carlo數值模擬結果．圖2顯示隨著α增大，曲線的峰值也逐漸增大，峰值的位置逐漸向左移．穩(wěn)態(tài)概率密度曲線的峰值越高，表示峰值處的概率越大，即分數階的階數α越大，系統(tǒng)在小振幅處取值的概率就越大．由此可以得出增大α可以減弱系統(tǒng)的穩(wěn)態(tài)響應．這個結論也可由圖3得到證實．由圖3可看出α增大時，位移的時間歷程圖變窄，反應系統(tǒng)的穩(wěn)態(tài)響應隨著α增大而減弱．另一方面我們也可以由廣義速度和位移的穩(wěn)態(tài)聯(lián)合概率密度來分析α對系統(tǒng)穩(wěn)態(tài)響應的影響．圖4給出了不同分數階α的穩(wěn)態(tài)聯(lián)合概率密度．可以看到隨著α的增大，圖像由粗矮狀逐漸變細、增高．這表示增大α時，系統(tǒng)的響應逐漸減弱．這里指出圖3和圖4中的參數取值和圖2中的參數取值是一樣的．

圖2 對于不同分數階數α幅值A的穩(wěn)態(tài)概率密度曲線Fig.2 Stationary probability density of amplitude A with different fractional order α

圖3（a） 對于不同分數階數α位移X的時間歷程圖Fig.3（a） Time history of the displacement X with different fractional order α

圖3（b） 對于不同分數階數α位移X的時間歷程圖Fig.3（b） Time history of the displacement X with different fractional order α

圖3（c） 對于不同分數階數α位移X的時間歷程圖Fig.3（c） Time history of the displacement X with different fractional order α

圖5給出對于不同分數階系數χ的幅值穩(wěn)態(tài)概率密度曲線．參數取值如下:α＝0．9，ω＝0．8，μ＝4，D＝0．04，δ1＝0．05，δ2＝ -0．05．其中實線表示理論結果，圓點表示原系統(tǒng)的Monte Carlo數值模擬結果．觀察圖5可以看出隨著χ的增大，概率密度曲線的峰值增大并且峰值的位置逐漸左移，由此分析得到增大χ同樣可以減弱系統(tǒng)的穩(wěn)態(tài)響應．

圖4（a） 關于廣義位移X和速度Y的穩(wěn)態(tài)聯(lián)合概率密度Fig.4（a） Joint stationary probability density of the generalized displacement X and velocity Y

圖4（b） 關于廣義位移X和速度Y的穩(wěn)態(tài)聯(lián)合概率密度Fig.4（b） Joint stationary probability density of the generalized displacement X and velocity Y

圖4（c） 關于廣義位移X和速度Y的穩(wěn)態(tài)聯(lián)合概率密度Fig.4（c） Joint stationary probability density of the generalized displacement X and velocity Y

在實際系統(tǒng)中噪聲對系統(tǒng)的影響是不可忽視的．圖6給出了不同噪聲強度D影響下幅值的穩(wěn)態(tài)概率密度曲線．實線表示理論結果，圓點表示原系統(tǒng)的Monte Carlo數值模擬結果．參數取值如下:α＝0．6，ω＝2，μ＝1．2，χ＝2，δ1＝0．05，δ2＝-0．05．觀察圖6可以看出增大噪聲強度D幅值A的穩(wěn)態(tài)概率密度曲線的峰值逐漸減小，并且峰值的位置逐漸向右移．說明增大D能夠增強系統(tǒng)的穩(wěn)態(tài)響應．

圖5 對于不同分數階系數χ幅值A的穩(wěn)態(tài)概率密度曲線Fig.5 Stationary probability density of the amplitude A with different fractional derivative intensitiesχ

圖6 對于不同噪聲強度D幅值A的穩(wěn)態(tài)概率密度曲線Fig.6 Stationary probability density of the amplitude A with different random intensities D

最后圖7給出了不同自然頻率ω影響下的穩(wěn)態(tài)概率密度．參數取值如下:α＝0．6，μ＝1．5，χ＝1．2，D＝0．06，δ1＝0．05，δ2＝-0．05．實線表示理論結果，圓點表示原系統(tǒng)的Monte Carlo數值模擬結果．從圖7中分析得增大ω可以減弱系統(tǒng)的穩(wěn)態(tài)響應．

圖7 對于不同自然頻率ω幅值A的穩(wěn)態(tài)概率密度曲線Fig.7 Stationary probability density of the amplitude A with different natural frequency ω

4 結論

本文研究了高斯白噪聲激勵下分數階Duffing-Van der Pol系統(tǒng)的穩(wěn)態(tài)響應.應用廣義諧和函數和隨機平均法將原系統(tǒng)轉化為關于幅值的平均伊藤微分方程，進而建立相應的PFK方程，并求解得到幅值的穩(wěn)態(tài)概率密度的近似解析解.對原系統(tǒng)進行Monte Carlo數值模擬驗證了理論結果的有效性.分析結果得出:增大分數階階數α、分數階系數χ以及自然頻率ω都能明顯地減弱系統(tǒng)的穩(wěn)態(tài)響應．相反地，增大高斯白噪聲的強度D可以增強系統(tǒng)的穩(wěn)態(tài)響應.

1 Wu W F，Lin Y K.Cumulant-neglect closure for non-linear oscillators under random parametric and external excitations.InternationalJournalofNon-LinearMechanics，1984，19（4）:349～362

2 Zeng Y，Zhu W Q.Stochastic averaging of n-dimensional non-linear dynamical systems subject to non-Gaussian wideband random excitations.International Journal of Non-Linear Mechanics， 2010，45（5）:572～586

3 Zhu W Q，Cai G Q.Random vibration of viscoelastic system under broad-band excitations.International Journal of Non-Linear Mechanics， 2011，46（5）:720～726

4 Ling Q，Jin X L，Huang Z L.Response and stability of SDOF viscoelastic system under wideband noise excitations.Journal of the Franklin Institute， 2011，348（8）:2026 ～2043

5 李東武，徐超.基于時頻域交替法的遲滯非線性振動系統(tǒng)的穩(wěn)態(tài)響應分析.動力學與控制學報，2016，14（3）:217～222 （Li D W， Xu C.Alternating time／frequency domain method for calculating the steady-state response of hysteresis nonlinear vibration systems.Journal of Dynamics and Control， 2016，14（3）:217～222 （in Chinese））

6 Mainardi F.Fractional calculus and waves in linear viscoelasticity:an introduction to mathematical models.Imperial College Press，2010:368

7 Alvelid M，Enelund M.Modelling of constrained thin rubber layer with emphasis on damping.Journal of Sound and Vibration， 2007，300（3-5）:662～675

8 Bagley R L， Torvik P J.Fractional calculus—A different approach to the analysis of viscoelastically damped structures.Aiaa Journal， 2012，21（5）:741～748

9 Bagley R L，Torvik P J.Fractional calculus in the transient analysis of viscoelastically damped structures.Aiaa Journal， 1985，23（6）:918～925

10 Machado J A T.And I say to myself:What a fractional world！.Fractional Calculus and Applied Analysis， 2011，14（4）:635～654

11 Machado J A T.Fractional order modelling of fractionalorder holds.Nonlinear Dynamics， 2012，70（1）:789～796

12 Machado J A T.Fractional Calculus:Application in modeling and control.Springer New York， 2013:279～295

13 Machado J A T，Costa A C，Quelhas M D.Fractional dynamics in DNA.Communications in Nonlinear Science and Numerical Simulation， 2011，16（8）:2963～2969

14 Mainardi F.Fractional relaxation-oscillation and fractional diffusion-wave phenomena.Chaos Solitons and Fractals，1996，7（9）:1461～1477

15 Chen L C，Wang W H，Li Z S，et al.Stationary response of Duffing oscillator with hardening stiffness and fractional derivative.International Journal of Non-Linear Mechanics，2013，48:44～50

16 Yang Y G，Xu W，Sun Y H，et al.Stochastic bifurcations in the nonlinear vibroimpact system with fractional derivative under random excitation.Communications in Nonlinear Science and Numerical Simulation， 2017，42:62～72

17 Chen L C，Zhu W Q.Stochastic jump and bifurcation of Duffing oscillator with fractional derivative damping under combined harmonic and white noise excitations.International Journal of Non-Linear Mechanics， 2011，46（10）:1324～1329

18 Huang Z L，Jin X L.Response and stability of a SDOF strongly nonlinear stochastic system with light damping modeled by a fractional derivative.Journal of Sound Vibration， 2009，319（3）:1121～1135

19 Yang Y G，Xu W，Gu X D，et al.Stochastic response of a class of self-excited systems with Caputo-type fractional derivative driven by Gaussian white noise.Chaos，Solitons and Fractals， 2015，77:190～204

20 Chen L C，Zhu W Q.First passage failure of SDOF nonlinear oscillator with lightly fractional derivative damping under real noise excitations.Probabilistic Engineering Mechanics， 2011，26（2）:208～214

21 Yang Y G，Xu W，Yang G D，et al.Response analysis of a class of quasi-linear systems with fractional derivative excited by Poisson white noise.Chaos， 2016，26（8）:083102

22 Xu Y，Li Y G，Liu D，et al.Responses of Duffing oscillator with fractional damping and random phase.Nonlinear Dynamics， 2013，74（3）:745～753

23 Xu Y，Li Y G，Liu D.Response of fractional oscillators with visc-oelastic term under random excitation.Journal ofComputational and Nonlinear Dynamics， 2014，9（3）:1081～1089

24 Spanos P D，Zeldin B A.Random vibration of systems with frequ-ency-dependent parameters or fractional derivatives.Journal of Engineering Mechanics， 1997，123（3）:290～292

25 Agrawal O P.Analytical solution for stochastic response of a fractionally damped beam.Journal of Vibration and A-coustics， 2004，126（4）:561～566

26 Ghorbanian P，Ramakrishnan S，Whitman A，et al.Nonlinear dynamic analysis of EEG using a stochastic Duffingvan der Pol oscillator model.ASME 2014 Dynamic Systems and Control Conference.American Society of Mechanical Engineers，2014:V002T16A001-V002T16A001

27 Xu Z，Cheung Y K.Averaging method using generalized harmonic functions for strongly non-linear oscillators.Journal of Sound Vibration， 1994，174（4）:563～576

28 Eugene W，Moshe Z.On the relation between ordinary and stochastic differential equations.International Journal of Engineering Science， 1965，3（2）:213～229

29 Diethelm K，Ford N J，Freed A D，et al.Algorithms for the fractional calculus:A selection of numerical methods.Computer Methods in Applied Mechanics and Engineering，2005，194（6-8）:743～773

30 Kai D， Ford N J， Freed A D.A predictor-corrector approach for the numerical solution of fractional differential equations.Nonlinear Dynamics， 2002，29（1-4）:3～22

STATIONARY RESPONSE OF DUFFING-VAN DER POL OSCILLATOR WITH FRACTIONAL DERIVATIVE UNDER GAUSSIAN WHITE NOISE?

Ma Yanyan Ning Lijuan?
（School of Mathematics and Information Science， Shaanxi Normal University， Xi′an 710119， China）

This paper investigates the stationary response of Duffing-Van der Pol oscillator with fractional derivative damping term under Gaussian white noise excitation.The corresponding Fokker-Planck-Kolmogorov（FPK）equation is firstly deduced by utilizing the stochastic averaging method and Stratonovich-Khasminskii theorem.And we then solve the FPK equation to obtain the stationary probability densities of system amplitude，which in fact can be used to describe the system response.It is then found that the reducing fractional derivative order enhances the system response， while the increasing fractional coefficient weakens the system response.Therefore，the fractional derivative damping term shows a great effect on the response of Duffing-Van der Pol oscillator.In addition， the response is also be influenced by other system parameters.Finally， the above analytical results reasonably fit with that of the Monte Carlo simulation.

response， fractional derivative， Duffing-van der Pol， Gaussian white noise， stochastic averaging method， Monte Carlo simulation

9 September 2016，revised 25 November 2016.

10.6052／1672-6553-2017-060

2016-09-09收到第1稿，2016-11-25收到修改稿.

?國家自然科學基金資助項目（11202120）和（GK201502007）

?通訊作者 E-mail:ninglijuan＠snnu.edu.cn

?The project supported by the National Natural Science Foundation of China（11202120） and （GK201502007）

?Corresponding author E-mail:ninglijuan＠snnu.edu.cn