ABSTRACT

This paper mainly analyzes the nonlinear random stability and random Hopf bifurcation of a four-dimensional high-speed rotor system under a random parameter excitation. The study of rotor system dynamics has made great progress in theory and practice. The system is affected by internal factors and external random wind and replaced by Gaussian color noise. By using the principle of random average, the quasi-hamiltonian system is convergent to a one-dimensional random ITO diffusion process, and then the maximum lyapunov exponential method is used to judge the local stability of the system, and the conditions for the local stability of the system are obtained. Then the Hopf bifurcation is simulated by the solution of FPK equation, namely the stationary probability density.

Keywords:Rotor System, Hopf Bifurcation, Stochastic Dynamics, Gauss Color Noise, Principle of Random Mean

一类转子系统的非线性随机稳定性及随机Hopf分岔

邓生文

兰州交通大学数理学院，甘肃 兰州

收稿日期：2020年3月25日；录用日期：2020年4月10日；发布日期：2020年4月17日

摘 要

本文主要分析一个随机参数激励下四维高速转子系统的非线性随机稳定性及随机Hopf分岔。转子系统动力学的研究在理论和实际操作也有了很大的进步。将系统受到的内部因素与外部随机风力影响用高斯色噪声代替。运用随机平均原理，将拟哈密顿系统收敛于一个一维伊藤随机扩散过程，然后运用最大李雅普诺夫指数法，来判断系统的局部稳定性，得到系统局部稳定的条件。然后通过FPK方程之解，即平稳概率密度来模拟系统发生Hopf分岔。

关键词 :转子系统，Hopf分岔，随机动力学，Gauss色噪声，随机平均原理

Copyright © 2020 by author(s) and Hans Publishers Inc.

This work is licensed under the Creative Commons Attribution International License (CC BY 4.0).

http://creativecommons.org/licenses/by/4.0/

1. 引言

转子–轴承系统是旋转机械中最关键的部分，而在设计和运转过程中，能够保证大型旋转机械中转子的稳定性和实用性是至关重要的，随着旋转机械向高效率、高标准发展，研究了不对称双圆盘转子–轴承系统在油膜力下的非线性耦合动力学行为。分析了非解析径向椭圆转子–轴承系统的稳定性和分岔，在稳定性分析中，运用打靶法和轨迹预测跟踪算法研究了系统线性不平衡响应。文献考虑了参数随机分布特性的转子动力学与可靠性 [1] [2] [3]。本节主要分析一个随机参数激励下四维高速转子系统的非线性随机稳定性及随机Hopf分岔。

2. 模型的建立

考虑一个带有广义非线性阻尼的转子模型 [4] [5]。

$\left(\begin{array}{cc}1& 0\\ 0& 1\end{array}\right)\left(\begin{array}{l}\ddot{x}\\ \ddot{y}\end{array}\right)+\chi \left(\begin{array}{cc}{\left|\stackrel{\dot{}}{x}\right|}^a+{\left|x\right|}^a& 0\\ 0& {\left|\stackrel{\dot{}}{y}\right|}^b+{\left|y\right|}^b\end{array}\right)\left(\begin{array}{l}\stackrel{\dot{}}{x}\\ \stackrel{\dot{}}{y}\end{array}\right)+\left(\begin{array}{cc}{k}_1/m& 0\\ 0& k_2/m\end{array}\right)\left(\begin{array}{l}x\\ y\end{array}\right)=(00)$

其中 $\chi$ 表示非线性阻尼系数， ${\omega }_n$ 是转子系统的固有频率。 ${\left|\stackrel{\dot{}}{x}\right|}^a+{\left|x\right|}^a$ 和 ${\left|\stackrel{\dot{}}{y}\right|}^b+{\left|y\right|}^b$ 表示非线性阻尼器的耗散效应，若 $a=b=0$，则系统的状态方程变成一个阻尼系数为 $2\chi$ 的形式，如果 $\chi$ 很小的情况下方程的近似解可以采用多种摄动方法得到，由于转子的质量相当大，所以 $\chi$ 的值可以被看作很小。

令 $a=b=2$，所以方程可以表示为：

$\{\begin{array}{l}\ddot{x}+\chi \left({\stackrel{\dot{}}{x}}^2+x^2\right)\stackrel{\dot{}}{x}+\left(k_1/m\right)x={\eta }_1\left(t\right)\\ \ddot{y}+\chi \left({\stackrel{\dot{}}{y}}^2+y^2\right)\stackrel{\dot{}}{y}+\left(k_1/m\right)y={\eta }_2(t)\end{array}$

考虑到外部的随机因素及其内部材质的不均质性，将系统收到的所有随机因素用Gauss色噪声代替，即 ${\eta }_i\left(t\right)$，$i=1,2$ 其统计性质满足以下条件：

$\langle {\eta }_i\left(t\right)\rangle =0$，$\langle {\eta }_i\left(t\right){\eta }_i\left(s\right)\rangle =\frac{D}{\tau }\exp \left(-\frac{\left|t-s\right|}{\tau }\right)$，$i=1,2$

其中，D和 $\tau$ 分别表示噪声强度和关联时间 [6]。

由于色噪声激励下的模型不是Markov过程，因此该模型的随机动力学行为的分析比较困难。本节将利用统一色噪声近似原理将模型简化成一个便于讨论的Gauss白噪声激励下的等效非线性模型。

2.1. 下面将色噪声模型进行等效白化 [6]

根据统一色噪声近似原理，色噪声激励下的模型可以由如下的白噪声激励下的方程组近似表示：

$\{\begin{array}{l}\ddot{x}+\chi \left({\stackrel{\dot{}}{x}}^2+x^2\right)\stackrel{\dot{}}{x}+\left(\frac{k_1}{m}\right)x={\eta }_1\left(t\right)\\ \ddot{y}+\chi \left({\stackrel{\dot{}}{y}}^2+y^2\right)\stackrel{\dot{}}{y}+\left(\frac{k_2}{m}\right)y={\eta }_2\left(t\right)\\ {\stackrel{\dot{}}{\eta }}_i\left(t\right)=-\frac{1}{\tau }{\eta }_i\left(t\right)+\frac{1}{\tau }{\xi }_i(t)\end{array}$

其中， $\xi \left(t\right)$ 是Gauss白噪声，满足 $\langle {\xi }_i\left(t\right)\rangle =0$，$\langle {\xi }_i\left(t\right){\xi }_i\left(s\right)\rangle =2D\delta \left(t-s\right)$，$i=1,2$。

为了后续计算方便，这里需要将方程组表示成一阶微分方程组的形式。故令

$x=q_1,\stackrel{\dot{}}{x}=p_1,y=q_2,\stackrel{\dot{}}{y}=p_2,{\omega }_1=k_1/m,{\omega }_2=k_2/m,\alpha =\chi$，可以得到

$\{\begin{array}{l}{\stackrel{\dot{}}{q}}_1=p_1\\ {\stackrel{\dot{}}{p}}_1=-{\omega }_1{q}_1-\alpha \left(p_1^2+q_1^2\right)p_1+{\eta }_1\left(t\right)\\ {\stackrel{\dot{}}{q}}_2=p_2\\ {\stackrel{\dot{}}{p}}_2=-{\omega }_2{q}_2-\alpha \left(p_2^2+q_2^2\right)p_2+{\eta }_2\left(t\right)\\ {\stackrel{\dot{}}{\eta }}_i\left(t\right)=-\left(1/\tau \right){\eta }_i\left(t\right)+\left(1/\tau \right){\xi }_i(t)\end{array}$

其中， $q_1,q_2$ 表示广义的位移， $p_1,p_2$ 表示广义的动量， $i=1,2$。

令 $f_1=-{\omega }_1{q}_1-\alpha \left(p_1^2+q_1^2\right)p_1$，$g_1=1$，$f_2=-{\omega }_2{q}_2-\alpha \left(p_2^2+q_2^2\right)p_2$，$g_2=1$

则方程组中(2)(4)式可以表示为

${\stackrel{\dot{}}{p}}_i=f_i+g_i\eta \left(t\right)$，$i=1,2$

将上式关于t求导，消去 $\eta \left(t\right)$ 可得：

${\ddot{p}}_i+c\left(p_i,\tau \right)\frac{{\stackrel{\dot{}}{p}}_i}{\tau }-\frac{f_i}{\tau }=\frac{g_i\xi \left(t\right)}{\tau }$

其中

$c_i\left(p_i,\tau \right)=1+\tau \frac{G_i}{g_i}$ , $G_i=f_i{g^{\prime }}_i-{f^{\prime }}_i{g}_i$ , $i=1,2$

应用统一色噪声近似方法，令 $p_i=0$，进行绝热消去，得到 ${\stackrel{\dot{}}{p}}_i={\bar{f}}_i+{\bar{g}}_i{\xi }_1\left(t\right)$，其中 ${\bar{f}}_i=\frac{f_i}{c}$，${\bar{g}}_i=\frac{g_i}{c}$，$i=1,2$ 故上式可以改写为

$\{\begin{array}{l}{\stackrel{\dot{}}{q}}_1=p_1\\ {\stackrel{\dot{}}{p}}_1=\frac{-{\omega }_1{\tau }^{-1}}{c_1\left(p_1,\tau \right)}{q}_1+\frac{-{\tau }^{-1}\alpha \left(p_1^2+q_1^2\right)}{c_1\left(p_1,\tau \right)}{p}_1+\frac{{\tau }^{-1}}{c_1\left(p_1,\tau \right)}{\xi }_1\left(t\right)\\ {\stackrel{\dot{}}{q}}_2=p_2\\ {\stackrel{\dot{}}{p}}_{12}=\frac{-{\omega }_2{\tau }^{-1}}{c_2\left(p_2,\tau \right)}{q}_1+\frac{-{\tau }^{-1}\alpha \left(p_2^2+q_2^2\right)}{c_2\left(p_2,\tau \right)}{p}_2+\frac{{\tau }^{-1}}{c_2\left(p_2,\tau \right)}{\xi }_2(t)\end{array}$

由于 $p_i$ 的非线性函数 $c_i\left(p_i,\tau \right)$，$i=1,2$ 位于模型第2、4个方程的分母上，致使白化模型是一个关于 $p_i$ 的强非线性模型，因而很难得到模型的最大Lyapunov指数解析表达式。此处考虑将模型部分线性化，得到一个等效非线性经济周期模型，从而降低模型的非线性强度和复杂性。

白化模型的等效非线性近似，可以得到模型的部分线性化等效近似表达式

$\{\begin{array}{l}{\stackrel{\dot{}}{q}}_1=p_1\\ {\stackrel{\dot{}}{p}}_1=-{\tau }^{-1}\left(k_{11}+k_{11}{p}_1\right)\left[{\omega }_1{q}_1+\alpha \left(p_1^2+q_1^2\right)p_1-{\xi }_1\left(t\right)\right]\\ {\stackrel{\dot{}}{q}}_2=p_2\\ {\stackrel{\dot{}}{p}}_2=-{\tau }^{-1}\left(k_{21}+k_{21}{p}_2\right)\left[{\omega }_2{q}_2+\alpha \left(p_2^2+q_2^2\right)p_2-{\xi }_2\left(t\right)\right]\end{array}$

记模型第2、4个方程与上述模型第2、4个方程的差为 $E_i$，则 $E_i$ 可以表示为

$E_i=-{\tau }^{-1}\left[\frac{1}{c_i\left(p_i,\tau \right)}-\left(k_{i1}+k_{i1}{p}_i\right)\right]\left({\omega }_i{q}_i+\alpha \left(p_i^2+q_i^2\right)p_i-2{\xi }_i\left(t\right)\right),i=1,2$

接着令 $E_i$ 的均值为零( ${\xi }_i\left(t\right)$ 值为零的白噪声)，得

$-{\tau }^{-1}\left[\frac{1}{c_i\left(p_i,\tau \right)}-\left(k_{i1}+k_{i1}{p}_i\right)\right]\left({\omega }_i{q}_i+\alpha \left(p_i^2+q_i^2\right)p_i\right)=0,i=1,2$

由于 ${\tau }^{-1}\left({\omega }_i{q}_i+\alpha \left(p_i^2+q_i^2\right)p_i\right),i=1,2$ 不可能恒为零，所以

整理上式，得

$\alpha k_{i2}{p}^4+\alpha k_{i1}{p}^3+\left({\tau }^{-1}{k}_{i2}+\omega k_{i2}+\alpha +\alpha q^2{k}_{i1}{k}_{i2}\right)p+{\tau }^{-1}{k}_{i1}+\omega k_{i1}+\alpha q^2{k}_{i1}-1=0$

考虑到在 $p=0,q=0$ 的邻域， $p^4$ , $p^3$ 与 $q^2$ 是零的高阶无穷小量，因此可以近似的认为上式等效于

$\left({\tau }^{-1}{k}_{i2}+\omega k_{i2}+\alpha \right)p+{\tau }^{-1}{k}_{i1}+\omega k_{i1}-1=0$

即

$\{\begin{array}{l}{\tau }^{-1}{k}_{i2}+\omega k_{i2}+\alpha =0\\ {\tau }^{-1}{k}_{i1}+\omega k_{i1}-1=0\end{array}$

求解上式，得

$\{\begin{array}{l}{k}_{i1}=\frac{1}{{\tau }^{-1}+\omega }\\ k_{i2}=-\frac{\alpha }{{\tau }^{-1}+\omega }\end{array}$

最后，代入原模型，得到该模型的等效非线性表达式为：

$\{\begin{array}{l}{\stackrel{\dot{}}{q}}_1=p_1\\ {\stackrel{\dot{}}{p}}_1=-a_{11}{q}_1+b_{11}\left(p_1^2+q_1^2\right)p_1+c_{11}{\xi }_1\left(t\right)\\ {\stackrel{\dot{}}{q}}_2=p_2\\ {\stackrel{\dot{}}{p}}_2=-a_{22}{q}_2+b_{22}\left(p_2^2+q_2^2\right)p_2+c_{22}{\xi }_2(t)\end{array}$

其中 $k_{i1}{\omega }_1=a_{ii}$，$k_{i2}\alpha =b_{ii}$，$k_{i1}+k_{i2}{p}_i=c_{ii}$，$i=1,2$。

该系统方程的哈密顿函数表示为：

$H=\frac{1}{2}\left(p_1{}^2+p_2{}^2\right)+\frac{1}{2}\left(a_{11}{q}_1^2+a_{22}{q}_2^2\right)$

根据拟不可积哈密顿函数 $H\left(t\right)$ 依概率1弱收敛与一个Itô随机扩散过程，则哈密顿函数的平均Itô随机微分方程 [8] 可以表示为：

$\text{d}\left(H\right)=\bar{m}\left(H\right)\text{d}t+\bar{\sigma }\left(H\right)\text{d}B(t)$

其中 $B\left(t\right)$ 是标准维纳过程, $m\left(H\right)$ 和 $\sigma \left(H\right)$ 是随机伊藤方程的漂移和扩散系数，根据随机平均原理：

由于高维系统在运用随机平均原理 [7] [8] 时需要多重积分，所以这里引入极坐标变换来计算漂移和扩散系数。

令 $q_1=\frac{R}{\sqrt{a_{11}}}\cos \theta ,q_{\text{2}}=\frac{R}{\sqrt{a_{22}}}\sin \theta$，

则方程转化为：

$\begin{array}{c}\bar{m}\left(H\right)=\frac{\text{2}\pi }{T\left(H\right)}{\displaystyle {\int }_{\Omega }[\left(b_{11}+b_{22}\right)A\left(H,\theta \right)-\left(\frac{b_{11}}{a_{11}}{\cos }^{\text{2}}\theta \text{+}\frac{b_{\text{22}}}{a_{22}}{\sin }^{\text{2}}\theta \right)B\left(H,\theta \right)}\\ +\frac{R^4}{2}\left(\frac{c_{11}^2\pi D_1}{a_{11}}{\cos }^2\theta +\frac{c_{22}^2\pi D_2}{a_{22}}{\sin }^2\theta \right)]\text{d}\theta \end{array}$

${\bar{\sigma }}^2\left(H\right)=\frac{4\pi }{T\left(H\right)}{\displaystyle {\int }_0^{\pi }\left(\frac{c_{11}^2\pi D_1}{a_{11}}{\cos }^2\theta +\frac{c_{22}^2\pi D_2}{a_{22}}{\sin }^2\theta \right)B\left(H,\theta \right)\text{d}\theta }$

$T\left(H\right)=\left(\frac{2\pi }{\sqrt{a_{11}{a}_{22}}}\right){\displaystyle {\int }_0^{\pi }{R}^2\text{d}\theta }$

$\Omega =\left\{\left(Q_1,\cdots ,Q_n,P_2,\cdots ,P_n\right)|H\left(Q_1,\cdots ,Q_n,0,P_2,\cdots ,P\right)\le H\right\}$

式中

$A\left(H,\theta \right)=H{R}^2-\frac{R^4}{4}\left(1+\frac{a_{11}+a_{22}}{2\sqrt{a_{11}{a}_{22}}}\sin 2\theta \right)-\frac{b_{11}+b_{22}}{24}{R}^6{\left(\frac{\cos \theta }{\sqrt{a_{11}}}-\frac{\sin \theta }{\sqrt{a_{22}}}\right)}^4$

$B\left(H,\theta \right)=\frac{H{R}^2}{2}-\frac{R^6}{6}\left(1+\frac{a_{11}+a_{22}}{2\sqrt{a_{11}{a}_{22}}}\sin 2\theta \right)-\frac{b_{11}+b_{22}}{32}{R}^8{\left(\frac{\cos \theta }{\sqrt{a_{11}}}-\frac{\sin \theta }{\sqrt{a_{22}}}\right)}^4$

而 $R^2$ 是下列方程之解：

$H-\frac{R^2}{2}\left(1+\frac{a_{11}+a_{22}}{2\sqrt{a_{11}{a}_{22}}}\sin 2\theta \right)-\frac{b_{11}+b_{22}}{8}{R}^4{\left(\frac{\cos \theta }{\sqrt{a_{11}}}-\frac{\sin \theta }{\sqrt{a_{22}}}\right)}^4=0$

2.2. 接下来主要讨论随机稳定性 [8]

由于 $\delta >1$，故当 $H=0$ 时，平均Itô随机微分方程的线性化形式为

$\text{d}H={\mu }_{1}H\text{d}t+{\mu }_2^{1/2}H\text{d}B(t)$

其中

${\mu }_1=-\frac{1}{2}\left(b_{11}+b_{22}\right)+\frac{1}{2}\left(\frac{c_{11}^2{D}_1}{a_{11}}+\frac{c_{22}^2{D}_2}{a_{22}}\right)\zeta$

${\mu }_2=\frac{1}{3}\left(\frac{c_{11}^2{D}_1}{a_{11}}+\frac{c_{22}^2{D}_2}{a_{22}}\right)\zeta$

$\zeta ={\displaystyle {\int }_0^{\pi }{\left(1+\frac{a_{11}+a_{22}}{2\sqrt{a_{11}{a}_{22}}}\sin \theta \right)}^{-2}\text{d}\theta }/{\displaystyle {\int }_0^{\pi }{\left(1+\frac{a_{11}+a_{22}}{2\sqrt{a_{11}{a}_{22}}}\sin \theta \right)}^{-1}\text{d}\theta }$

接着可以得到，当 $H\to 0$ 时， $\bar{m}\left(H\right)$ 和 ${\bar{\sigma }}^2\left(H\right)$ 的渐进形式

$\bar{m}\left(H\right)={\mu }_{1}H+o\left(H\right)$，${\bar{\sigma }}^2\left(H\right)={\mu }_2{H}^2+o(H2)$

基于Oseledec乘积遍历性定理 [8]，可以得到模型的最大Lyapunov指数为：

$\begin{array}{l}\lambda =\underset{t\to \infty }{\lim }\frac{1}{t}\ln H^{\frac{1}{2}}\left(t\right)\\ =\frac{1}{2}\left[m^{\prime }\left(0\right)-{\left({\sigma }^2(0)\right)}^2\right]\\ =-\frac{1}{4}\left(b_{11}+b_{22}\right)+\frac{1}{6}\left(\frac{c_{11}^2{D}_1}{a_{11}}+\frac{c_{22}^2{D}_2}{a_{22}}\right)\zeta \end{array}$

系统平凡解保持局部渐近稳定的充要条件为： $\lambda <0$，即

$\frac{2}{3}\left(\frac{c_{11}^2{D}_1}{a_{11}}+\frac{c_{22}^2{D}_2}{a_{22}}\right)\zeta <\left(b_{11}+b_{22}\right)$

3. 随机分岔分析

由Itô随机微分方程解得相应的FPK方程为：

$\frac{\partial p}{\partial t}=-\frac{\partial }{\partial H}\left[a\left(H\right)p\right]+\frac{1}{2}\frac{{\partial }^2}{\partial H^2}\left[b\left(H\right)p\right]$

式中

$a\left(H\right)=\bar{m}\left(H\right),b\left(H\right)={\bar{\sigma }}^2(H)$

因为当 $H\to 0$ 时， $\bar{m}\left(H\right)=0$，${\bar{\sigma }}^2\left(H\right)=0$，所以根据边界类别原理可得，左边界 $H\to 0$ 属于第一类奇异边界，进而可以得到漂移指数 ${\alpha }_H$ 和扩散指数 ${\beta }_H$，以及特征标值 $c_H$

${\alpha }_H=1$，${\beta }_H=2$，

$c_H=\underset{x\to x_H^{+}}{\lim }\frac{2a\left(x\right){\left(x-x_H\right)}^{{\alpha }_H-{\beta }_H}}{b\left(x\right)}=-\frac{3}{2\zeta }\left(b_{11}+b_{22}\right)+\frac{3}{2}\left(\frac{c_{11}^2{D}_1}{a_{11}}+\frac{c_{22}^2{D}_2}{a_{22}}\right)$

则平稳概率密度函数为：

$f\left(H\right)=C{H}^{v_1}\exp \left[u_{1}H\right]$

$v_1=c_H-{\alpha }_H=-2+\left(-\frac{3}{2\zeta }\left(b_{11}+b_{22}\right)+\frac{3}{2}\left(\frac{c_{11}^2{D}_1}{a_{11}}+\frac{c_{22}^2{D}_2}{a_{22}}\right)\right)$

$u_1=\frac{3}{2}\left(\frac{c_{11}^2{D}_1}{a_{11}}+\frac{c_{22}^2{D}_2}{a_{22}}\right)$

因此可得联合概率密度函数为：

$\begin{array}{l}f\left(q_1,q_2,p_1,p_2\right)=C{\left(\frac{1}{2}\left(p_1^2+p_2^2\right)+\frac{1}{2}\left(a_{11}{q}_1^2+a_{22}{q}_2^2\right)\right)}^{v_1}\\ \times \exp \left[u_1\left(\frac{1}{2}\left(p_1^2+p_2^2\right)+\frac{1}{2}\left(a_{11}{q}_1^2+a_{22}{q}_2^2\right)\right)\right]\end{array}$

4. 数值模拟

下面主要分析随机Hopf分岔，通过模拟平稳概率密度函数 $f\left(H\right)$ 和联合概率密度函数 $f\left(q_1,q_2,p_1,p_2\right)$ 的图像随着参数变化而变化。选定 $b_{11}+b_{22}$ 作为分岔参数，然后令 $C=1$，$a_{11}=0.2$，$a_{22}=0.3$，$c_{11}=c_{22}=1$，$D_1=D_2=-1$，$\zeta =1$。为了能够模拟联合概率密度函数，同样令 $q_2=p_2=0.2$。

图1. 当 $b_{11}+b_{22}=-9.67$ 时， $v_1=0.01$ 系统平稳概率密度函数和联合概率密度函数

图2. 当 $b_{11}+b_{22}=-9.89$ 时， $v_1=0.335$ 系统平稳概率密度函数和联合概率密度函数

图3. 当 $b_{11}+b_{22}=-11.67$ 时， $v_1=3$ 系统平稳概率密度函数和联合概率密度函数

平稳概率密度函数与联合概率密度函数图像如图1~3所示，当参数 $b_{11}+b_{22}$ 在−9.89附近变化时，平稳概率密度图像从单调递减变为单峰状函数，系统发生Hopf分岔；而联合概率密度函数则随着的逐渐变大由单峰状函数变为火山口状，系统发生Hopf分岔。

5. 小结

在本文中，分析了一个非线性阻尼转子系统的随机稳定性及Hopf分岔。首先在一个四维转子–轴承系统中考虑了随机高斯色噪声，然后将高斯色噪声进行白化，再运用哈密顿系统理论对模型进行稳定性分析。在随机平均的过程中引入了极坐标变换，得到了伊藤随机微分方程的漂移系数和扩散系数，通过最大李雅普诺夫指数法，计算出该系统的稳定性条件；最后研究了系统随机Hopf分岔行为。在数值时，只考虑其中2个变量，另外2个变量令为参数的形式，得到联合概率密度函数和平稳概率密度函数随着分分岔参数变化的图像，得以验证得到的理论。

文章引用

邓生文. 一类转子系统的非线性随机稳定性及随机Hopf分岔
Nonlinear Stochastic Stability and Stochastic Hopf Bifurcation of a Class of Rotor Systems[J]. 应用数学进展, 2020, 09(04): 501-508. https://doi.org/10.12677/AAM.2020.94061

参考文献