周期與色噪聲聯合作用下分數階Duffing振子非平穩響應的無記憶方法

摘要 基于統計線性化提出了一種求解周期與色噪聲激勵聯合作用下分數階Duffing系統非平穩響應的無記憶方法。將系統響應分解為確定性周期和零均值隨機分量之和，則原非線性運動方程可等效地化為一組耦合的、分別以確定性和隨機動力響應為未知量的分數階微分方程。利用無記憶化方法將確定性和隨機分數階微分方程轉化為相應的常微分方程。利用統計線性化方法處理隨機常微分方程，得到關于隨機響應二階矩的李雅普諾夫方程。利用數值算法聯立求解李雅普諾夫微分方程和確定性常微分方程。通過Monte Carlo模擬，驗證此方法的適用性和精度。

關鍵詞 非平穩響應; 分數階系統; 無記憶方法; 統計線性化; 聯合激勵

引 言

分數階微積分近幾十年在工程界得到廣泛的關注［1］。分數階微積分的一個重要工程應用是黏彈性材料的力學模型的建立。與標準線性固體模型（Standard Linear Solid Model）相比，分數階導數模型能以較少的參數擬合實驗獲得黏彈性松馳和蠕變數據［2］；此外，分數階微分方程可很好地描述動力激勵下裝備黏彈性控制裝置（如：天然橡膠支座［3］、黏滯阻尼器［4］和黏彈性阻尼器［5?6］等）的動力行為。因此，研究這類分數階運動方程的求解方法成為學者們關注的重點。然而，大量研究集中在分數階運動方程的確定性分析［7?8］，關于它們的隨機分析尚較少涉及。分數階線性系統的隨機動力響應可采用頻域方法［9］和時域方法［10］，它們可視為整數階系統經典線性隨機方法的延伸。同樣地，分數階非線性系統的隨機分析方法，如統計線性化［11］、隨機平均法［12?13］、路徑積分［14］等也可視為整數階系統經典非線性方法的擴展。最近，Kong等在Monte Carlo模擬的框架下提出了分數階線性或非線性系統平穩響應的多諧波平衡方法［15］和非平穩響應的小波?Galerkin法［16］。此外，還有一些非經典方法，如利用特征向量展開求解具有1/2階分數階阻尼的單自由度系統的確定性［17］和隨機響應［18］。

與整數階導數不同，分數階導數的一個重要特性在于其記憶性，即當前狀態的分數階導數依賴于所有歷史狀態，它導致了求解確定性分數階微分方程的長持時響應非常耗時。為此，人們發展了一種處理分數階導數的無記憶法，即采用變量代換將分數階導數等效地轉化為若干整數階導數。這種方 法最早由Yuan等［19］在確定性分數階系統動力響應分析的背景下提出。隨后，大量研究集中于解釋該方法的等效力學模型［20］、提高計算精度［21?22］和改進計算效率［23］等方面。文獻［24］將這種方法應用于線性隨機振動；文獻［25］將這種方法應用于非線性系統隨機動力響應的Monte Carlo模擬。

另一方面，實際工程中存在系統或結構受到確定性和隨機激勵聯合作用的情況。研究者們主要關注白噪聲作用下非線性系統的平穩響應密度函數解［26］和系統穩定性［27］。最近，考慮到工程隨機動力系統的特殊性，人們發展了聯合激勵下整數階非線性集中質量系統平穩響應的統計線性化方法［28］，并將其擴展應用于非線性分數階集中質量［29］和連續質量系統［30］。這種方法本質上基于統計線性化，具有處理多種工程隨機問題的潛力。隨后，孔凡等將這種方法推廣到聯合激勵下多自由度滯回系統［31?32］、非平穩響應［33］的情況。

本文利用無記憶方法處理分數階導數，并基于統計線性化發展一種聯合激勵下分數階Duffing振子非平穩響應的方法。方法的具體步驟如下：首先，將響應分解為周期和隨機響應分量的組合，得到兩個耦合的、關于周期和隨機響應分量的分數階微分方程；隨后，利用無記憶方法將分數階微分方程去記憶化后轉化為不含分數階項的確定性和隨機常微分方程組；接著，利用統計線性化方法處理隨機常微分方程組，得到關于隨機響應二階矩的李雅普諾夫方程；最后，同時求解確定性常微分方程組和李雅普諾夫方程可得所有響應分量。

4 數值算例

不失一般性，確定性激勵采用單諧波形式，即：

式中 Fs為確定性激勵幅值，ω0為確定性激勵頻率。隨機激勵為均勻調制白噪聲或均勻調制色噪聲。

4.1 Qs（t）為白噪聲

4.1.1 典型參數設置的情況

選取系統參數m=1， c=0.4， k=1，α=0.5， ε=0.1，確定性激勵參數Fs=0.5，ω0=1，隨機激勵參數S0=0.4/π，A=1， μ1=0.1， μ2=0.2。其中m=1，k=1和S0=0.4/π由歸一化處理而來；Fs選為0.5，以使確定性響應與隨機響應在一個數量級；另外，為使系統處于共振情況，將確定性激勵頻率ω0選為1；拉蓋爾節點數n=2，n=4和n=8。利用獨立高斯分布序列生成白噪聲樣本10000條，并乘以調制函數得到均勻調制白噪聲。由所建議方法和蒙特卡羅模擬（Monte Carlo Simulation， MCS）得到位移均值與標準差，分別如圖1（a）和圖1（b）所示。

從圖1 （a）可以看出，即使當n=2時，本文所建議方法得到的響應均值與蒙特卡羅模擬所得結果吻合較好。由圖1（b）可知，建議方法所得的標準差精度隨著n值的增加而提升；當n=8時，此方法所得到的標準差與MCS結果完全吻合。此外，由圖1（b）可見，所建議方法可以很好地捕捉由于突加激勵、諧波激勵的耦合效應和調制函數引起的非平穩效應。

T過大會導致隨機分量標準差的時間平均過小且隨各參數變化不明顯，因此需對T的取值加以限制。選定標準差對時間積分占總積分面積97%左右的時間點為截止時間T。本文中50 s時標準差積分占1000 s時的97.34%，所以截止時間選為50 s。

4.1.2 非線性強度的影響

其他參數與4.1.1節選取相同，非線性強度系數ε取0～1。采用所建議方法與蒙特卡羅模擬計算得到的確定性響應模的時間平均對比如圖2所示；標準差時間平均對比如圖3所示。

從圖2， 3可見，在其他參數不變的情況下，確定性響應隨非線性強度增加而降低；位移標準差隨非線性強度增加而略有降低。在所有非線性強度情況下，所建議方法的精度隨著n的增加而提高；n=8時，兩種方法吻合較好。

4.1.3 激勵幅值的影響

為研究確定性激勵幅值對此方法適用性的影響，保持其他參數不變，確定性激勵幅值Fs取0～2。本文所建議方法與蒙特卡羅模擬計算得到的確定性響應模的時間平均對比如圖4所示；時變標準差的時間平均對比如圖5所示。

由圖4，5可知，在其他參數不變的情況下，確定性響應隨諧波激勵幅值增大而增大；位移標準差隨諧波激勵幅值增大而減小。不同確定性激勵幅值下，拉蓋爾節點數對隨機響應的影響比對確定性響應的影響更加明顯。n=8時，所建議方法得到的結果與蒙特卡羅模擬高度吻合，顯示了該方法在不同確定性響應幅值下的適用性。

4.1.4 確定性諧波激勵頻率的影響

為研究確定性激勵頻率對此方法適用性的影響，保持其他參數不變，確定性激勵頻率ω0取0.1～2。采用本文建議方法與蒙特卡羅模擬得到的確定性響應分量的模的時間平均對比如圖6所示；隨機響應非平穩標準差的時間平均如圖7所示。

由圖6，7可見，諧波激勵頻率對確定性響應和隨機響應均有較大影響，存在使他們達到極值的諧波激勵頻率。其中，使確定性響應和隨機響應達到最大的諧波激勵頻率約為1.2 rad/s。圖6類似純諧波激勵下非線性振子的“幅頻響應曲線”，對于本文考慮的硬化Duffing系統呈不對稱峰值形態。拉蓋爾節點數的增加會明顯提高所建議方法的精度；當n=8時，在共振和非共振頻率下本文建議方法與蒙特卡羅模擬所得結果均吻合較好。

4.1.5 隨機激勵強度的影響

為研究隨機激勵強度對此方法適用性的影響，保持其他參數不變，功率譜強度S0取10?4～102。采用本文所建議方法與蒙特卡羅模擬得到的確定性響應模的時間平均對比如圖8所示；隨機響應時變標準差的時間平均對比如圖9所示。

由圖8， 9可知，隨機激勵功率譜強度對確定性響應和隨機響應均有一定影響：隨機激勵譜強度的增大會使確定性響應減小，并使隨機響應增大。拉蓋爾節點數的增加會提高所建議方法的精度；當n=8時，所建議方法與蒙特卡羅模擬所得結果吻合較好。結合4.1.3節中得到的結果可知，確定性或隨機激勵幅值會分別正向影響確定性響應或隨機響應，而分別反向影響隨機響應或確定性響應。

綜上，對于分數階Duffing系統在確定性諧波與調制白噪聲聯合作用下的響應，本文所建議的方法在考察的不同參數設置情況下，均有良好的適用性。

4.2 Qs（t）為色噪聲

4.2.1 典型參數設置的情況

選取系統參數m=1，c=0.4，k=1，ε=0.1， α=0.5；確定性激勵參數Fs=0.4，ω0=1；隨機激勵參數A=1， μ1=0.1，μ2=0.2，ζg=0.5， ωg=1， S0=0.1。利用譜表現方法生成色噪聲樣本10000條，乘以調制函數后得到均勻調制色噪聲。建議方法和蒙特卡羅模擬得到位移均值與標準差分別如圖10（a）圖10（b）所示。

從圖10（a），（b）可見，在考慮的參數設置情況下，本文建議方法的精度隨n值的增大而提高；當n=8時，建議方法與蒙特卡羅模擬高度吻合。在隨機響應方面，所建議方法能很好地捕捉突加激勵、確定性激勵分量和隨機激勵分量調制函數導致的非平穩性。

下面分析其他參數對建議方法適用性的影響。確定性響應模的時間平均和隨機響應標準差時間平均的定義同4.1.1節。

4.2.2 非線性強度的影響

保持其他參數不變，非線性強度系數ε取0～1。采用建議方法與蒙特卡羅模擬計算得到的確定性響應模的時間平均對比如圖11所示；標準差時間平均對比如圖12所示。

由圖11，12可見，在其他參數不變的情況下，確定性和隨機響應均隨非線性強度增加而降低。本文建議方法的精度隨n值的增大而提高；對于確定性響應，當n=8時，所建議方法與MCS高度吻合；對于隨機響應，n=8時，所建議方法具有滿意的精度。以上結果表明了所建議方法在不同非線性強度下的適用性。

4.2.3 激勵幅值的影響

為研究確定性激勵幅值對此方法適用性的影響，保持其他參數不變，確定性激勵幅值Fs取0～2。本文所建議方法與蒙特卡羅模擬計算得到的確定性響應模的時間平均對比如圖13所示；時變位移標準差的時間平均值對比如圖14所示。

由圖13，14可見，在其他參數不變的情況下，確定性響應隨諧波激勵幅值增大而增大，而位移標準差隨諧波激勵幅值增大而減小。后者顯示了確定性分量與隨機分量之間的耦合性：確定性激勵大小會影響隨機響應量大小。此外，拉蓋爾節點數n的增加會提升所建議方法的精度；當n=8時，建議方法與蒙特卡羅模擬所得結果高度吻合。以上分析均反映了建議方法在不同確定性激勵幅值情況下的適用性。

4.2.4 確定性諧波激勵頻率的影響

為研究確定性激勵頻率對此方法適用性的影響，保持其他參數不變，確定性激勵頻率ω0取0.1～2。本文所建議方法與蒙特卡羅模擬得到的確定性響應模的時間平均對比如圖15所示；非平穩隨機響應標準差的時間平均如圖16所示。

由圖15，16可見，諧波激勵頻率對確定性和隨機響應有較大影響，存在使它們達到極值的諧波激勵頻率。確定性響應與隨機響應達到極值的頻率約為1.2 rad/s。圖15所示的曲線與非線性系統在諧波激勵下的“幅頻響應曲線”類似，呈現“硬特性”系統特有的峰值不對稱特征。圖16顯示了隨機響應受到確定性激勵頻率的影響，表明了二者之間是耦合的。拉蓋爾節點數n的增加會提升所建議方法的精度；當n=8時，所建議方法與蒙特卡羅模擬所得結果高度吻合。以上分析均表明了建議方法在不同確定性激勵頻率下的適用性。

4.2.5 隨機激勵強度的影響

同樣，為研究色噪聲激勵強度對此方法適用性的影響，保持其他參數不變，功率譜強度S0取10?4～102。采用本文所建議方法與蒙特卡羅模擬得到的確定性響應模的時間平均對比如圖17所示；隨機響應時變標準差的時間平均對比如圖18所示。

由圖17，18可見，隨機激勵功率譜強度對確定性和隨機響應均有一定影響：隨著隨機激勵功率譜強度的增大，確定性響應逐漸減小，而隨機響應逐漸增大。結合4.2.3節中的結果可知，確定性或隨機激勵幅值會分別正向影響確定性或隨機響應幅值，反向影響隨機或確定性響應幅值；因此，確定性和隨機分量之間是相互耦合的。此外，建議方法的精度隨著n值的增大而提高；對于確定性響應而言，n=8時，建議方法與蒙特卡羅模擬所得結果高度吻合。以上分析均反映了所建議方法在不同隨機激勵強度下的適用性。

綜上所述，對于分數階Duffing系統在確定性諧波與隨機激勵（白噪聲或色噪聲）聯合作用下的響應，本文所建議方法在不同參數設置的情況下，均有良好的適用性。所建議方法的計算效率較蒙特卡羅模擬有顯著優勢。以調制白噪聲與確定性激勵聯合作用下的響應計算為例，使用所建議方法計算拉蓋爾節點數n=8時用時僅需0.39 s，而10000個樣本的蒙特卡羅模擬需要39.24 s。由此可見，所建議方法在保證了較高精確性的前提下，計算效率上也顯著提高。

5 結論與展望

對單自由度分數階Duffing振子在周期和非平穩隨機激勵聯合作用下的響應進行了研究，得到的結論如下：

（1） 基于無記憶化方法提出了一種求解周期和非平穩隨機激勵聯合作用下單自由度分數階Duffing振子非平穩響應的統計線性化方法，并利用蒙特卡羅法驗證了該方法的適用性與精度；

（2） 算例表明，所提方法在考慮的參數設置范圍內具有較好的精度，同時其計算精度隨拉蓋爾節點數的增加而增加；

（3） 在多數情況下，當拉蓋爾節點數為8時，即可獲得較為理想的結果，甚至在共振頻率下也具有較為理想的精度；

（4） 所提方法在保證較高精度的前提下，計算效率相較于蒙特卡羅法得到了顯著的提升。

參考文獻

1Cai M， Li C. Numerical approaches to fractional integrals and derivatives： a review［J］. Mathematics， 2020， 8（1）： 43.

2Di Paola M， Pirrotta A， Valenza A. Visco-elastic behavior through fractional calculus： an easier method for best fitting experimental results［J］. Mechanics of Materials， 2011， 43（12）： 799-806.

3Koh C G， Kelly J M. Application of fractional derivatives to seismic analysis of base‐isolated models［J］. Journal of Earthquake Engineering and Structural Dynamics， 1990， 19（2）： 229-241.

4Makris N， Constantinou M C. Fractional-derivative Maxwell model for viscous dampers［J］. Journal of Structural Engineering， 1991， 117（9）： 2708-2724.

5Zhao X， Wang S， Du D， et al. Optimal design of viscoelastic dampers in frame structures considering soil-structure interaction effect［J］. Shock and Vibration， 2017， 2017： 9629083.

6Lewandowski R， Pawlak Z. Dynamic analysis of frames with viscoelastic dampers modelled by rheological models with fractional derivatives［J］. Journal of Sound and Vibration， 2011， 330（5）： 923-936.

7孔凡， 侯召旭， 徐軍， 等. 基于多諧波平衡法的滯回分數階系統穩態動力響應［J］. 振動工程學報， 2021， 34（3）： 552-558.

Kong F， Hou Z， Xu J， et al. Steady-state response determination of a hysteretic system endowed with fractional elements via a multi-harmonic balance method［J］. Journal of Vibration Engineering， 2021， 34（3）： 552-558.

8Di Paola M， Pinnola F P， Spanos P D. Analysis of multi-degree-of-freedom systems with fractional derivative elements of rational order［C］. ICFDA'14 International Conference on Fractional Differentiation and Its Applications 2014， IEEE， 2014：1-6.

9Zeldin B A， Spanos P D. Random vibration of systems with frequency-dependent parameters or fractional derivatives［J］. Journal of Engineering Mechanics， 1997， 123（3）： 290-292.

10Agrawal O P. Stochastic analysis of dynamic systems containing fractional derivatives［J］. Journal of Sound and Vibration， 2001， 247（5）： 927-938.

11Spanos P D， Evangelatos G I. Response of a non-linear system with restoring forces governed by fractional derivatives—time domain simulation and statistical linearization solution［J］. Soil Dynamics and Earthquake Engineering， 2010， 30（9）： 811-821.

12Chen L， Wang W， Li Z， et al. Stationary response of Duffing oscillator with hardening stiffness and fractional derivative［J］. International Journal of Non-Linear Mechanics， 2013， 48： 44-50.

13孫春艷， 徐偉. 含分數階導數項的隨機Duffing振子的穩態響應分析［J］. 振動工程學報， 2015， 28（3）： 374-380.

Sun C Y， Xu W. Stationary response analysis for a stochastic Duffing oscillator comprising fractional derivative element［J］. Journal of Vibration Engineering， 2015， 28（3）： 374-380.

14Matteo A D， Kougioumtzoglou I A， Pirrotta A， et al. Stochastic response determination of nonlinear oscillators with fractional derivatives elements via the Wiener path integral［J］. Probabilistic Engineering Mechanics， 2014， 38： 127-135.

15Kong F， Spanos P D. Response spectral density determination for nonlinear systems endowed with fractional derivatives and subject to colored noise［J］. Probabilistic Engineering Mechanics， 2020， 59： 103023.

16Kong F， Zhang Y， Zhang Y. Non-stationary response power spectrum determination of linear/non-linear systems endowed with fractional derivative elements via harmonic wavelet［J］. Mechanical Systems and Signal Processing， 2022， 162（168）： 108024.

17Suarez L E， Shokooh A. An eigenvector expansion method for the solution of motion containing fractional derivatives［J］. Journal of Applied Mechanics， 1997， 64（3）： 629-635.

18Agrawal O P. Stochastic analysis of a 1-D system with fractional damping of order 1/2［J］. Journal of Vibration and Acoustics， 2002， 124（3）： 454-460.

19Yuan L， Agrawal O P. A numerical scheme for dynamic systems containing fraction derivatives［J］. Journal of Vibration and Acoustics， 2002， 124（2）： 321-324.

20Schmidt A， Gaul L. On a critique of a numerical scheme for the calculation of fractionally damped dynamical systems［J］. Mechanics Research Communications， 2006， 33（1）： 99-107.

21Agrawal O P. A numerical scheme for initial compliance and creep response of a system［J］. Mechanics Research Communications， 2009， 36（4）： 444-451.

22Agrawal O P. A modified memory-free scheme and its Simulink implementation for FDEs［J］. Physica Scripta， 2009， 2009： 014031.

23Liu Q X， Chen Y M， Liu J K. An improved Yuan?Agrawal method with rapid convergence rate for fractional differential equations［J］. Computational Mechanics， 2019， 63（4）： 713-723.

24Di Paola M， Failla G， Pirrotta A. Stationary and non-stationary stochastic response of linear fractional viscoelastic systems［J］. Probabilistic Engineering Mechanics， 2012， 28： 85-90.

25Failla G， Pirrotta A. On the stochastic response of a fractionally-damped Duffing oscillator［J］. Communications in Nonlinear Science and Numerical Simulation， 2012， 17（12）： 5131-5142.

26Chen L C， Zhu W Q. Stochastic averaging of strongly nonlinear oscillators with small fractional derivative damping under combined harmonic and white noise excitations［J］. Nonlinear Dynamics， 2009， 56（3）： 231-241.

27Chen L C， Zhu W Q. Stochastic stability of Duffing oscillator with fractional derivative damping under combined harmonic and white noise parametric excitations ［J］. Acta Mechanica， 2009， 207（1）： 109-120.

28Spanos P D， Zhang Y， Kong F. Formulation of statistical linearization for MDOF systems subject to combined periodic and stochastic excitations［J］. Journal of Applied Mechanics， 2019， 86（10）： 101003.

29孔凡， 晁盼盼， 徐軍， 等. 隨機與諧和聯合激勵下分數階非線性系統的統計線性化方法［J］. 振動工程學報， 2021， 34（4）： 756?764.

Kong F， Chao P P， Xu J， et al. Statistical linearization method for nonlinear Duffing oscillator under combined random and harmonic excitations［J］. Journal of Vibration Engineering， 2021， 34（4）： 756?764.

30Spanos P D， Malara G. Nonlinear vibrations of beams and plates with fractional derivative elements subject to combined harmonic and random excitations［J］. Probabilistic Engineering Mechanics， 2020， 59： 103043.

31Kong F， Spanos P D. Stochastic response of hysteresis system under combined periodic and stochastic excitation via the statistical linearization method［J］. Journal of Applied Mechanics， 2021， 88（5）： 051008.

32孔凡， 韓仁杰， 張遠進， 等. 色噪聲與確定性諧波聯合激勵下Bouc-Wen動力系統響應的統計線性化方法［J］. 振動工程學報， 2022， 35（1）： 82?92.

Kong F， Han R J， Zhang Y J， et al. Stochastic response of a hysteresis system subjected to combined periodic and colored noise excitation via the statistical linearization method［J］. Journal of Vibration Engineering， 2022， 35（1）： 82?92.

33孔凡， 韓仁杰， 張遠進. 確定性周期與隨機激勵聯合作用下非線性系統非平穩響應的統計線性化方法［J］. 振動工程學報， 2022， 35（3）： 625-634.

Kong F， Han R J， Zhang Y J. Non-stationary response of non-linear systems subjected to combined periodic and non-stationary stochastic excitation via the statistical linearization method［J］. Journal of Vibration Engineering， 2022， 35（3）： 625?634.

34Roberts J B， Spanos P D. Random Vibration and Statistical Linearization［M］. Dover Publications， 2003.