孔凡 韩仁杰 张远进 李书进

摘要：提出了一种用于求解色噪声和确定性谐波联合作用下单自由度Bouc-Wen系统响应的统计线性化方法。基于系统响应可分解为确定性谐波和零均值随机分量之和的假定，将原滞回运动方程等效地化为两组耦合的且分别以确定性和随机动力响应为未知量的非线性微分方程。利用谐波平衡法求解确定性运动方程，利用统计线性化方法求解色噪声激励下的随机运动方程。由此，可导出关于确定性谐波响应分量F ourier级数和随机响应分量二阶矩的非线性代数方程组。利用牛顿迭代法对上述耦合的代数方程组进行求解。数值算例验证了此方法的适用性和精度。

关键词：统计线性化;Bouc-Wen滞回模型;谐波平衡法;联合激励;牛顿迭代法

中图分类号：0324;TU311.4

文献标志码：A

文章编号：10044523（ 2022）01-0082-11

DOI： 10.1638 5/j .cnki.issn.10044523.2022.01.009

引言

随机振动分析方法已被广泛地应用于工程科学的各个领域。由Booton[1]和Caughey[2]先后提出的统计线性化（ Statistical Linearization，SL）方法是解决非线性系统随机振动常用的方法之一[3]。该方法同样适用于分数阶非线性系统[4]。最近，基于小波分析时域一频域联合分辨的概念[5]，作者与其合作者提出了时一频域等效线性化方法，并将其应用于完全非平稳随机过程激励下的非线性系统[6]。关于统计线性化方法最新进展的综述，可参阅文献[7]。

然而，某些情况下，工程结构会同时受到确定性周期和随机激励作用。例如，旋转式飞机[8]经常受到色噪声和谐波激励联合作用;风力发电机的叶片对湍流的响应[9]等。因此，谐波与随机激励联合作用下非线性系统响应的研究越来越受到广大学者的关注[10-11]。在此背景下，人们提出了几种解析和数值方法。这些方法通常利用各种确定性方法与随机方法的组合求解耦合的确定性与随机微分方程。其中，包括确定性线性化和高斯线性化或矩截断方法的组合[12-14]、多尺度法与高斯线性化或矩截断方法的组合[15-16]、多尺度法与随机平均法的组合[17]、谐波平衡法与随机平均法的组合[l8-19]、确定性平均法和统计线性化或高斯矩截断方法的组合[20]谐波平衡法与高斯线性化或矩截断方法的组合[21-22]、随机平均法与统计线性化的组合[23]。此外，还可利用基于马氏随机过程的方法求解响应概率密度函数，以及考察联合激励下非线性系统的跳跃、分岔现象。即通过数值方法（如中心差分法[24]和路径积分法[25]）或解析方法[26]求解随机平均法得到的FP（ Fokker-Planck）方程或CK（ Chapman-Kolmogorov）方程;抑或直接根据原随机动力系统的Ito随机微分方程，用路径积分法[27]或胞映射法[28]求解响应的概率密度函数。

从前面的文献综述可以看出，几乎所有研究者都关注多项式非线性系统。例如Duffing[21，25]，Vander Pol[14.18]和Duffing-Rayleigh[15]振子。然而，非线性多项式并不能準确地描述材料在大变形情况下的滞回现象，即材料或构件的本构关系或力一位移曲线依赖于它的加载历程。就本文作者所知，极少有研究者关注滞回系统在随机与谐和联合激励作用下的响应。然而，在很多工程实际中却会出现这种情况，如近断层地震作用下的铅芯橡胶隔震结构。

本文提出一种求解色噪声和确定性谐波联合作用下单白由度B ouc-Wen系统响应的统计线性化方法。该方法基于系统响应可分解为确定性谐波和零均值随机分量之和的假定。基于该假定，可将原滞回运动方程等效地化为两组耦合、分别以确定性和随机动力响应为未知量的非线性微分方程。随后，利用谐波平衡法求解确定性运动方程，并利用统计线性化方法求解色噪声激励下的随机运动方程。由此，可导出关于确定性谐波响应分量Fourier级数和随机响应分量二阶矩的非线性代数方程组。利用牛顿迭代法对上述耦合的代数方程组进行求解。最后，数值算例验证此方法的适用性。

1 动力学方程

单自由度Bouc-Wen系统在确定性谐波和随机色噪声联合激励下的运动方程为：

因此，谐波和随机联合激励下的原运动方程（式（1））和滞回方程（式（2》可转化为确定性微分方程（式（9）和（12））和随机微分方程（式（10）和（13）），且二者之间是耦合的。下节中，将利用谐波平衡法求解确定性分量的Fourier系数。

2 谐波响应分量的谐波平衡法

结合式（16），（17）和（27），（28）可求解确定性响应Fourier级数Co，Do，Uo，Vo。然而，上述方程中除未知响应Fourier系数外还耦合有未知随机响应特征值（ρ，σx，σz）。因此，还需要更多代数方程使上述方程组完备。下节中，将对式（IO）和（13）使用统计线性化方法以得到ρ，σx，σz;与C0，D0，U0，V0之间的其他代数关系。

3 随机响应分量的统计线性化方法

将式（34）-（37）代入式（32）和（33），可知等效线性参数是随时间呈谐和变化的，因为其中含有均值过程μx和μz（见式（14）和（7））。考虑到当t→∞。时标准差σx与σz趋于循环平稳（cyclo-stationary），可消除由μx和μz引起的快变性。所以，等效线性化参数可取一个周期（T0=2π/ω0）内的平均值，即：

从上述分析可见，等效线性参数ce和ke由7个未知量C0，D0，U0，V0，σx，σz和ρ确定。可通过随机振动的状态空间法得出随机参数σx，σz，ρ等效线性化参数ce和ke的联系。

作为演示，假定随机激励的功率谱密度为[30]：

本文利用软化和硬化B ouc-Wen系统验证所提出方法的适用性。软化B ouc-Wen系统的滞回参数取A=1，γ= 0.5，β= 0.5，n=1'a=0.1;硬化Bouc-Wen系统β= -0.35，），γ=0.65。本文所提方法与Monte Carlo模拟（Monte Carlo simulation，MCS）的结果对比如图1，2所示。其中，MCS中，样本激励由谱表现方法生成。图1是软化B ouc-Wen系统响应的对比结果，可见本文提出方法得到的响应均值和标准差与10000个样本的MCS所得结果总体吻合良好。注意到，MCS的方差在达到平稳后仍出现类似简谐的抖动，这是由于确定性和随机响应耦合效应造成的，其中还包含响应样本统计的随机性因素。对达到平稳后的MCS方差进行若干整数周期上的时间平均可得到响应方差呈现谐波变化的基线值，以下误差分析均以该基线值为标准。图1（a）所示的位移响应对比中，二者所得均值幅值相差 2.65%;图l（b）中，方差达到平稳后呈谐和变化的幅值较小，所建议方法得到的平稳方差与MCS估计的平稳方差相差约-13.85%。本文提出的方法涉及对等效线性参数的时间平均，抹去了响应方差的谐波变化特征，值得进一步改进。试算表明：谐波频率一定时，幅值越大，或谐波幅值一定時，频率越接近共振频率，响应方差的谐波变化特征越明显。

同样地，所建议方法对硬化系统也有很好的计算精度，如图2所示。具体而言，位移响应均值的幅值相差- 0.43%，平稳方差相差-12.19%。以上误差均在一般统计线性化方法的合理误差范围之内。

5.1简谐激励频率的影响

需要注意的是，在推导式（18），（19）和（34）～（37）的过程中，假定了 的值为小量。当

的值由0到1逐渐增大时，式（18），（19）和（34）-（37）的近似值与其精确值之间相差如附录中图Al--A3所示。当谐波激励频率接近系统白振频率时，或谐波激励幅值增大时，谐波响应分量

会增大。因此，讨论谐波激励在不同幅值与频率下方法的适用性是非常重要的。图3和4为谐波激励幅值F0=0.8时，响应标准差（

）和确定性响应幅值随简谐激励频率的变化曲线。可见，该情况下由本文所建议方法求得的响应标准差与10000个样本Monte Carlo模拟所得到的结果符合较好。

Monre Carlo模拟结果表明，当激励的频率等于非线性结构的白振频率时，

和

达到峰值。当F0= 0.8时，二者的峰值分别为2.23和1.20。其中，两种方法得出的响应最大差别如表l所示。由表可知，在整个频率范围内，F0=0.8时，所建议方法的最大误差均在一般统计线性化方法的合理误差范围内。

进一步研究此方法对于硬化Bouc-Wen系统的适用性。同样地，图5和6为谐波激励幅值为F0=0.3时，随机响应分量标准差（ ）和谐和响应分量幅值随谐波激励频率变化的曲线以及所建议方

可见，本文所建议方法与Monte Carlo模拟值在多数情况下吻合较好。当F0= 0.3时，和 的峰值分别为0.93和0.96。两种方法所得结果的最大误差列于表2中。结果表明，F0=0.3时的最大误差均在一般统计线性化方法误差的合理范围内。

5.2 简谐激励幅值的影响

显然，谐波激励的幅值影响

）和

）的大小，从而进一步影响所建议方法的精度。就此，采用ω0=1（共振），讨论所建议方法精度与简谐激励幅值的关系。对软化Bouc-Wen系统，图7和8分别为随机动力响应分量的标准差（

）和谐和响应分量幅值在简谐激励频率为ω0=1时，随简谐激励幅值变化的曲线。

计算表明，指标

）随简谐激励幅值单调变化。当ω0=1时，二指标最大值分别为2.52和1.05。图7，8同时给出了共振频率下，指标小于阈值1时，简谐激励幅值的范围。当二指标均小于预定阈值时，可视为满足本文所设假定条件。此时，将所建议方法得到的结果与MonteCarlo模拟之间最大相对误差列于表3中。结果表明，满足本文所设假定条件时，建议方法的误差均在一般统计线性化方法误差的合理范围内。此外，简谐激励幅值等于0时，对应系统处于完全随机激励的情况。由图7可知，ω0=1时，系统随机位移分量的标准差随着简谐激励幅值增大而增大，随机速度和滞回位移分量的标准差随简谐激励幅值增大而减小。由图8可知，ω0=1时确定性总位移、速度和滞回位移幅值随简谐激励幅值增大而增大。当简谐激励幅值处于假定应用范围时，本文所建议方法和MC模拟得到的确定性响应幅值之间的差别极小;简谐激励幅值不处于假定应用范围时，本文所建议方法也能准确地捕捉上述趋势。

同样地，对于硬化Bouc-Wen系统，图9和10分别为随机响应分量标准差（ ）和确定性响应幅值随简谐激励幅值的变化曲线。图中均给出了两种方法在激励频率为ω0=1下的响应对比。

指标

的值随简谐激励幅值增大而增大，相应地，所建议方法的精度变差。当ω0=1时，指标最大值分别为3.50和3.43。同样地，图9和10给出了指标满足预定阈值时的简谐激励幅值区间。此时，将两种方法所得结果的最大相对误差列于表4中。指标满足本文所做假定时，建议方法与Monte Carlo对比的相对误差在一般统计线性化方法的合理误差范围内。同样地，随机响应分量标准差和确定性响应分量幅值均随简谐响应幅值有各白的变化趋势，本文所建议方法均能在假定适用范围内较好地捕捉这一趋势。

6 结论

本文提出了一种求解Bouc-Wen滞回系统在确定性谐波与色噪声联合激励作用下的统计线性化方法。该方法基于系统响应可分解为确定性谐波和零均值随机分量之和的假定。基于该假定，将原滞回运动方程等效地化为了以确定性和随机动力响应为未知量的两组耦合的非线性微分方程。随后，利用谐波平衡法求解了确定性运动方程，并利用统计线性化方法求解了色噪声激励下的随机运动方程。由此，导出了关于确定性谐波响应分量Fourier级数和随机响应分量二阶矩的非线性代数方程组。利用牛顿迭代法求解了上述耦合的代数方程组。最后，数值算例验证了此方法的适用性。考察了软化B ouc-Wen系统和硬化Bouc-Wen系统在不同激励幅值和共振与非共振情况下的响应。结果表明了几乎在所有满足适用性条件的情况下，此方法都有合理的精度。注意到，本文提出采用统计线性方法求解联合激励下滞回系统的随机动力响应，与基于马尔可夫过程的方法相比，在牺牲了一定精度的情况下，大大拓展了该方法的适用性范围。因此，更适合于求解工程随机动力系统近似响应。

参考文献：

[1] Booton R C.The analysis of nonlinear control systemswith random inputs[C].Proceedings of the Symposiumon Nonlinear Circuit Analysis. New York. 1953：369-391.

[2] Caughey T K.Equivalent linearization techniques[Jl.Journal of the Acoustical Society of America. 1962. 35（11）：1706-1711.

[3]

Spanos P D.Stochastic linearization in structural dynamics[J].Applied Mechanics Reviews， 1981. 34（1）： 1-8.

[4]孙春艳，徐伟.含分数阶导数项的随机Duffing振子的稳态响应分析[J].振动T程学报，2015，28（3）：42- 48.

Sun Chunyan， Xu Wei. Stationary response analysis fora stochastic Duffing oscillator comprising fractional derivative element [J] . Journal of Vibration Engineering，2015.28（ 3） ：4248.

[5]孔凡 ，王恒 ，徐军 ，等 .分数阶导数系统响应功率谱密度 的小波-Galerkin方法 [J].振动工.程学报 ， 2018， 31（ 4） ： 671-680.

Kong Fan， Wang Heng， Xu Jun. et al.Wavlet-Galerkin based approach for determining powerspectral density of a system endowed with fractional derivative damping [ J] . Journal of Vibration Engineering ，2018 . 31（ 4） ： 671-680.

[6]Kong F， Kougioumtzoglou I A， Spanos P D. et al.Nonlinear system response evolutionary power spectraldensity determination via a harmonic wavelets based（Jalerkin technique [J] . International Journal for Multiscale Computational Engineering， 2016. 14 （3） ：255272.

[7]Elishakoff Isaac， Crandall Stephen H. Sixty years ofstochastic linearization technique [ J] . Meccanica ， 2017 .52（ 1-2） ： 299305.

[8]Hatchell B K. Mauss F J， Amaya I A， et al. Missilecaptive carry monitoring and helicopter identification us-ing a capacitive microelectromechanical systems accelerometer [Jl. Structural Health Monitoring， 2012， 11（ 2） ： 213224.

[9]Megerle B， Stephen Rice T. McBean I. et al. Numerical and experimental investigation of the aerodynamicexcitation of a model lowpressure steam turbine stageoperating under low volume flow [Jl. Journal of Engineering for （Jas Turbines and Power， 2013. 135（1） ：012602.

[10]Harne R L， Dai Q. Characterizing the robustness andsusceptibility of steady-state dynamics in postbuckledstructures to stochastic perturbations [J]. Journal ofSound and Vibration. 2017 . 395（ 2） ： 258271.

[11]陳林聪，朱位秋，谐和与宽带噪声联合激励下含分数导数型阻尼 的 Duffing振子 的平稳 响应 [J] .力学学报 ， 2010 ， 27（ 3） ： 517-521.

Chen Lincong ， Zhu Weiqiu. Stationary response of duff-ing oscillator with fractional derivative damping undercombined harmonic and wide band noise excitations [J] .Chinese Journal of Applied Mechanics， 2010， 27 （ 3） ：517-521.

[12]Ellermann K. On the determination of nonlinear re-sponse distributions for oscillators with combined harmonic and random excitation [Jl. Nonlinear Dynamics，2005. 43（ 3） ： 305318.

[13]Budgor A B. Studies in nonlinear stochastic processes.iii. approximate solutions of nonlinear stochastic differential equations excited by CJaussian noise and harmonicdisturbances [ J] . Journal of Statistical Physics ， 1977 ， 17（ 1） ： 21-44.

[14]Bulsara A R， Lindenberg K， Shuler K E. Spectral anal-ysis of a nonlinear oscillator driven by random and periodic forces. I. Linearized theory [Jl. Journal of Statistical Physics ， 1982 ， 27（ 4） ： 787-808.

[15]Nayfeh A H. Serhan S J. Response statistics ofnon-linear systems to combined deterministic and ran-dom excitations [ J] . Pergamon ， 1990 ， 25（ 5） ： 493-509.

[16]戎海武，王向東，孟光，等，简谐与随机噪声联合激励下Van der P01-Duffing系统的响应 [Jl.振动工程学报 ，2003 ， 16（ 4） ： 502-505.

Rong Haiwu， Wang Xiangdong， Meng Guang， et al.Response of Van der P01-Duffing oscillator to combinedexcitation of harmonic and random inputs [ J] . Joumal ofVibration Engineering ， 2003 ， 16（ 4） ： 502505.

[17]Rong H， Wang X， Xu W. et al. Resonant response of anon-linear vibroimpact system to combined determinis-tic harmonic and random excitations [Jl. InternationalJournal of Non-Linear Mechanics， 2010， 45 （ 5） ：474481.

[18]Manohar C S， Iyengar R N. Entrainment in van derPol's oscillator in the presence of noise [ J] . InternationalJournal

of Non-Linear Mechanics ，

1991，

26 （ 5） ：679686.

[19]Haiwu R， Wei X， （Juang M， et al. Response of a Duffing oscillator to combined deterministic harmonic andrandom excitation [ J ] . Journal of Sound and Vibration .2001， 242（ 2） ： 362368.

[20]Anh N D ， Hieu N N. The Duffing oscillator under com-bined periodic and random excitations [J] . ProbabilisticEngineering Mechanics ， 2012 ， 30（ 2） ： 27-36.

[21]Zhu H. Guo S. Periodic response of a Duffing oscillatorunder combined harmonic and random excitations [Jl.Journal of Vibration and Acoustics， Transactions of theASME. 2015. 137（4） ： 041015.

[22]Spanos P D ， Zhang Y ， Kong F. Formulation of statistical linearization for MD-（）F systems subject to com-bined periodic and stochastic excitations [Jl. Journal ofApplied Mechanics ， 2019 ， 86（ 10） ： 101003.

[23]Anh N D， Zakovorotny V L， Hao D N. Response anal-ysis of van der Pol oscillator subjected to harmonic andrandom excitations [Jl. Probabilistic Engineering Mechanics， 2014， 37（5）： 51-59.

[24]Narayanan S . Kumar P. Numerical solutions of FokkerPlanck equation of nonlinear systems subjected to ran-dom and harmonic excitations [Jl. Probabilistic Engineering Mechanics ， 2011 . 27（ 1） ： 3546.

[25]Huang Z L， Zhu W Q， Suzuki Y. Stochastic averagingof strongly nonlinear oscillators under combined harmon-ic and white-noise excitations [Jl. Journal of Sound andVibration ， 2000 ， 238（ 2） ： 233256.

[26]Cai G Q， Lin Y K. Nonlinearly damped systems undersimultaneous broadband and harmonic excitations [Jl.Nonlinear Dynamics ， 1994， 6（ 2） ： 163-177.

[27]Yu J S， Lin Y K. Numerical path integration of a non-homogeneous Markov process [ J] . International Journalof Non-Linear Mechanics ， 2004 . 39（ 9） ： 1493-1500.

[28]Han Q， Xu W， Sun J. Stochastic response and bifurcation of periodically driven nonlinear oscillators by thegeneralized cell mapping method [ J] . Physica A ： Statistical Mechanics and Its Applications， 2016， 458 （4） ：115-125.

[29]Wen Y K. Method for random vibration of hystereticsystems [ J] . Joumal of the Engineering Mechanics Division， 1976 ， 102（ 2） ： 249263.

[30]Kanai K. Semiempirical formula for the seismic charac-teristics of the ground [ J] . Transactions of the Architec-tural Institute of Japan . 1957 . 35.

[31]Jameson A. Solution of the equation AX+XB= C by in-version of an MXM or NXN matrix [ J] . SIAM Joumalon Applied Mathematics ， 1968 . 16（ 5） ： 1020-1023.