DSC  >> Vol. 9 No. 1 (January 2020)

    一类环上的三自由度耦合van der Pol方程的同伦分析方法
    The Homotopy Analysis Method for a Class of Three-Degree-of-Freedom Coupled van der Pol Oscillators

  • 全文下载: PDF(2436KB) HTML   XML   PP.22-39   DOI: 10.12677/DSC.2020.91003  
  • 下载量: 113  浏览量: 150   国家自然科学基金支持

作者:  

胡佩倩,钱有华:浙江师范大学,数学与计算机科学学院,浙江 金华

关键词:
同伦分析方法耦合van der Pol振子周期解Homotopy Analysis Method Coupled van der Pol Oscillator Periodic Solution

摘要:

本文主要运用同伦分析方法研究一类三自由度耦合van der Pol振子环系统周期解的近似表达式。我们先给出此系统的一般方程,再将此一般方程分成四类来讨论。第一,所有振子都同步运动;第二,两个振子同步运动,而第三个振子以一无关的方式运动(除它与第二个振子有相同周期的振动外);第三,环上相邻的振子之间彼此都相位差1/3周期的运动;第四,两个振子相位差1/2周期,而第三个振子2倍于它们的频率振动。利用四种不同类型的van der Pol振子来说明同伦分析方法的有效性与广泛应用性,且将此方法与数值积分法进行了比较,结果发现得到的解析解与数值解具有很高的吻合性。

In this paper, the approximate expression of periodic solution for a class of 3-DOF coupled van der Pol oscillator ring system is studied by homotopy analysis method. We first give the general equation of the system, and then divide the general equation into four categories: First, all oscillators move synchronously; second, two oscillators move synchronously, while the third vibrator moves in an independent manner (except that it oscillates in the same period as the second); third, the phase difference between adjacent oscillators on the ring is one third-period motion. Fourth, the two oscillators are 1/2 cycle apart, and the third vibrates at twice their frequency. Four different types of van der Pol oscillators are used to illustrate the validity and wide application of the homotopy analysis method.

1. 引言

随着科学技术的发展,越来越多实际问题显现出非线性的本质,故而非线性问题早已渗透于我们身边的方方面面。对于非线性问题,我们通常转变成非线性方程来解决。而对于非线性方程的求解,不像线性问题可以用线性叠加原理来求解。非线性问题求解起来相对比较复杂,而对于此问题,国内外专家学者近年来也进行了诸多研究。对于强非线性系统的定量分析方法主要有L-P摄动法 [1] 、同伦分析方法 [2] [3] [4] [5] 、双曲函数摄动法 [6] 等。

廖世俊教授于1998年首次提出同伦分析方法,并将其应用到非线性振动的研究中 [7] [8],文献 [9] 给出非线性单自由度振动系统的两个周期公式,该公式对系统任何的振幅都是一致有效的。廖世俊等学者 [10] 应用同伦分析方法,举例并分析研究了非牛顿磁流体流动系统的解析解。钱有华等人 [11] 对同伦分析方法进行了改造,得到一类强非线性耦合van der Pol-Duffing系统的解析近似周期解。2016年,崔继峰 [12] 利用同伦分析方法并结合多尺度思想,对受迫van der Pol-Duffing振子的周期解等性质进行了研究。2016年,许天亮等人 [13] 利用同伦分析方法研究了几类简单非线性微分方程问题,并验证了同伦分析方法在求解部分非线性微分方程的优越性。2017年,郭江明 [14] 用同伦分析方法探究了带有时滞的耦合van der Pol系统双Hopf分岔和周期解。2018年,杨孟 [15] 等学者采用同伦分析方法研究了圆盘上非牛顿幂律流体的Marangoni边界层流动传热问题,利用求解化简后的非线性问题,从而获得一系列的近似解析解。

耦合van der Pol方程也是社会诸多学者研究的一个热点问题,其在动力学性质方面取得许多显著的成果。2007年,Rompala Kevin等 [16] 基于3-DOF耦合van der Pol系统,研究系统的动态行为。2009年,Jiashi Tang等 [17] 研究了对耦合van der Pol极限环的振幅进行控制的有效方法,并分析了部分弱非线性系统,得出了控制参数同振幅之间的关系,为如何有效控制振幅提出有效建议。同年,钱有华等 [18] 巧妙利用同伦分析方法探究了参数激励的非自制薄板系统的周期解。2016年,吴俊 [19] 利用KAM理论分析证明拟合周期解存在于2-DOF耦合van der Pol方程中。2018年,石兆羽等 [20] 利用混沌振子对微弱信号敏感以及对强噪声具有良好免疫力的特性,提出基于耦合van der Pol-Duffing振子系统检测微弱信号的新方法。2019年,施添添等 [21] 通过讨论特征方程根分布情况研究含时滞的大规模van der Pol-Duffing耦合振子系统的同步和异步周期振荡、概周期运动以及混沌吸引子等非线性动力学现象。

本文研究的耦合van der Pol振子环结构如图1所示 [22] [23]。而环状结构在医学、光学、电路、机电传感器和神经网络等诸多领域有着重要的应用。2007年,Fuad,Baolong Lu和Surendra Singh [24] 研究了带调制泵的双向固态环形激光器的准同步、同步和非同步混沌振荡的机制,并用最大李亚普诺夫指数进行了计算。2010年,Fernando,Ana Paula和Carla [25] 研究了由两个(单向或双向)环形细胞通过缓冲细胞耦合而成的耦合细胞网络的动力学特性。2016年,Faralli [26] 等人论证了在具有母线和环形拓扑结构的光子集成片上网络中,在同一微环中添加或删除两个波长相同的光信号的可行性。2019年,Xu,Yu和Wang [27] 为了避免简单环拓扑结构的缺点,提出了一种具有时延的耦合双环网络。采用分治算法分析相邻矩阵特征值的分布研究了周期动力学平衡和层间同步的稳定性,并进一步分析了耦合环与反射波之间的动力学相互作用。

本文主要建立三自由度耦合van der Pol振子环的精确近似解析解,且利用此系统的周期解证明了同伦分析方法的有效性和巨大的应用潜力。此外,还对同伦分析方法和数值积分法的计算结果进行了比较。结果表明,在时间历程图中,即使时间t的值在很大的区间内变化,同伦分析方法的解析解与数值积分解也能很好地吻合。

Figure 1. The ring of three-degree-of-freedom van der Pol oscillators with coupled

图1. 耦合三自由度van der Pol振子环

2. 同伦分析方法简述

我们先引入多自由度耦合振子的一般表达式:

M q ¨ + G q ˙ + K q = F ( q ) , (2-1)

其中q为n维未知向量,即 q = ( x 1 , x 2 , , x n ) T ,且 q ¨ 为q关于时间变量t的二阶微分,类似地, q ˙ 为q关于时间变量t的一阶微分,M,G,K为n阶方阵,而耦合项 F ( q ) 则是关于变量 q ˙ , q , t 的向量值函数。

首先,我们先引入一个非线性算子:

N [ u ( r , t ) ] = M 2 u ( r , t ) t 2 + G u ( r , t ) t + K u ( r , t ) F [ u ( r , t ) ] , (2-2)

其中

u ( r , t ) = ( x 1 ( t ) , , x n ( t ) ) T , (2-3)

u ( r , t ) t = ( d x 1 d t , , d x n d t ) T , (2-4)

2 u ( r , t ) t 2 = ( d 2 x 1 d t 2 , , d 2 x n d t 2 ) T , (2-5)

而这里的 u ( r , t ) 是未知向量函数,r和t分别是空间与时间变量。

由同伦分析方法可知,我们需要构造一零阶形变方程

( 1 p ) { L [ Φ ( r , t , p ) u 0 ( r , t ) ] } = p h H ( t ) N [ Φ ( r , t , p ) ] , (2-6)

其中 p [ 0 , 1 ] 为嵌入变量,L为辅助线性算子, u 0 ( r , t ) 为初始猜测解,而h和 H ( t ) 分别为辅助参数和辅助函数。辅助参数h对级数解的收敛性有决定性作用,故我们可画出h-ω曲线来确定h的值来确保级数解收敛。

易知,当 p = 0 时,零阶形变方程(2-6)演变为 Φ ( r , t , 0 ) = u 0 ( r , t ) ;而当 p = 1 时,(2-6)又可化简为 Φ ( r , t , 1 ) = u ( r , t ) 。所以我们可得出结论:当嵌入变量p由 0 1 的变化过程中,解 Φ ( r , t , p ) 从初始猜测解 u 0 ( r , t ) 变化到精确解 u ( r , t )

u m ( r , t ) = 1 m ! m Φ ( r , t , p ) p m | p = 0 (2-7)

我们把 Φ ( r , t , p ) 展开成关于嵌入变量p的Taylor级数,可得

Φ ( r , t , p ) = u 0 ( r , t ) + m = 1 + u m ( r , t ) p m . (2-8)

若选取适当的辅助线性算子L,辅助参数h和辅助函数 H ( t ) ,再结合有效的初始猜测解便可使得幂级数(2-8)在 p = 1 时收敛,因而可得

u ( r , t ) = u 0 ( r , t ) + m = 1 + u m ( r , t ) . (2-9)

向量 u m 定义为

u m = { u 0 ( r , t ) , u 1 ( r , t ) , , u m ( r , t ) } . (2-10)

我们对零阶形变方程(2-6)关于嵌入变量p求m次导数,再将求导得到的式子左右两边同除以 m ! ,并令 p = 0 ,则我们就可得到m阶形变方程

L [ u m ( r , t ) χ m u m 1 ( r , t ) ] = h H ( t ) R m ( u m 1 , r , t ) , (2-11)

其中

χ m = { 0 , m 1 1 , m > 1 (2-12)

R m ( u m 1 , r , t ) = 1 ( m 1 ) ! m 1 N ( Φ ( r , t , p ) ) p m 1 | p = 0 . (2-13)

原则上m阶形变方程是线性的,而对于线性问题我们可借助数学软件如Mathematica求解。

3. 同伦分析方法的应用

本论文所研究的系统是三自由度耦合van der Pol方程,其一般表达式为

μ 1 x ¨ 1 + ε 1 ( x 1 2 λ ) x ˙ 1 + η 1 x 1 = f 1 ( x 1 , x 2 , x 3 ) μ 2 x ¨ 2 + ε 2 ( x 2 2 λ ) x ˙ 2 + η 2 x 2 = f 2 ( x 1 , x 2 , x 3 ) μ 3 x ¨ 3 + ε 3 ( x 3 2 λ ) x ˙ 3 + η 3 x 3 = f 3 ( x 1 , x 2 , x 3 ) (3-1)

其中 f 1 , f 2 , f 3 互为耦合项。对于上述所给的一般方程,我们不妨规定 μ 1 = μ 2 = μ 3 ε 1 = ε 2 = ε 3 η 1 = η 2 = η 3 ,而对于方程(3-1)中 f i ( i = 1 , 2 , 3 ) 的取值,我们将根据振子运动模式分成四类来讨论 [28] :第一,所有振子都同步运动;第二,两个振子同步运动,而第三个振子以一无关的方式运动(除它与第二个振子有相同周期的振动外);第三,环上相邻的振子之间彼此都相位差1/3周期的运动;第四,两个振子相位差1/2周期,而第三个振子2倍于它们的频率振动。我们将针对上述的四种振子不同的运动状态研究系统的周期解问题。在文章最后,我们会采用数值模拟的形式对其进行验证,表明改进的同伦分析方法作用于此类系统是可取的,并希望该方法在其他科学领域能得到广泛应用。

要解决上述问题,我们先考虑如下的三自由度非线性耦合van der Pol振子:

μ 1 x ¨ 1 + ε 1 ( x 1 2 λ ) x ˙ 1 + η 1 x 1 = f 1 ( x 1 , x 2 , x 3 ) μ 2 x ¨ 2 + ε 2 ( x 2 2 λ ) x ˙ 2 + η 2 x 2 = f 2 ( x 1 , x 2 , x 3 ) μ 3 x ¨ 3 + ε 3 ( x 3 2 λ ) x ˙ 3 + η 3 x 3 = f 3 ( x 1 , x 2 , x 3 ) (3-2)

其中 x ¨ i ( i = 1 , 2 , 3 ) 为关于变量t的二阶微分, x ˙ i ( i = 1 , 2 , 3 ) 为关于变量t的一阶微分, x i ( i = 1 , 2 , 3 ) 为关于变量t的未知实函数, f i ( i = 1 , 2 , 3 ) 是耦合项 μ i , ε i , η i , λ ( i = 1 , 2 , 3 ) 为已知实参数。

我们设该系统的周期解频率为 ω ,然后引入一个新变量 τ = ω t ,将该新变量代入系统(3-2),得到

μ 1 ω 2 u ¨ 1 + ε 1 ω ( u 1 2 λ ) u ˙ 1 + η 1 u 1 = f 1 ( u 1 , u 2 , u 3 ) μ 2 ω 2 u ¨ 2 + ε 2 ω ( u 2 2 λ ) u ˙ 2 + η 2 u 2 = f 2 ( u 1 , u 2 , u 3 ) μ 3 ω 2 u ¨ 3 + ε 3 ω ( u 3 2 λ ) u ˙ 3 + η 3 u 3 = f 3 ( u 1 , u 2 , u 3 ) (3-3)

方程满足初始条件

u 1 ( 0 ) = a , u 1 ( 0 ) = b , u 2 ( 0 ) = c , u 2 ( 0 ) = d , u 3 ( 0 ) = e , u 3 ( 0 ) = 0 (3-4)

其中导数表示对变量 τ 微分。

若系统(3-3)的周期解可以用基函数 { cos ( k τ ) , sin ( k τ ) | k = 0 , 1 , 2 , } 表示,则有

u 1 ( τ ) = k = 0 + ( a 1 , k cos k τ + b 1 , k sin k τ ) u 2 ( τ ) = k = 0 + ( a 2 , k cos k τ + b 2 , k sin k τ ) u 3 ( τ ) = k = 0 + ( a 3 , k cos k τ + b 3 , k sin k τ ) (3-5)

我们可以假定初始猜测解可表示为以下形式

u 1 , 0 ( τ ) = a 0 cos τ + b 0 sin τ , u 2 , 0 ( τ ) = c 0 cos τ + d 0 sin τ , u 3 , 0 ( τ ) = e 0 cos τ . (3-6)

之后,我们还需要引入辅助线性算子

L ( Φ 1 ( τ , p ) Φ 2 ( τ , p ) Φ 3 ( τ , p ) ) = ω 0 2 ( 2 Φ 1 ( τ , p ) τ 2 + Φ 1 ( τ , p ) 2 Φ 2 ( τ , p ) τ 2 + Φ 2 ( τ , p ) 2 Φ 3 ( τ , p ) τ 2 + Φ 3 ( τ , p ) ) , (3-7)

线性算子L满足下列性质

L ( C 1 sin τ + C 2 cos τ C 3 sin τ + C 4 cos τ C 5 sin τ + C 6 cos τ ) = 0. (3-8)

对于系统(3-2),我们再定义非线性算子

N ( Φ 1 ( τ , p ) Φ 2 ( τ , p ) Φ 3 ( τ , p ) ) = ( μ 1 Ω 2 ( p ) 2 Φ 1 ( τ , p ) τ 2 + ε 1 Ω ( p ) [ Φ 1 2 ( τ , p ) λ ] Φ 1 ( τ , p ) τ + η 1 Φ 1 ( τ , p ) f 1 ( Φ 1 ( τ , p ) , Φ 2 ( τ , p ) , Φ 3 ( τ , p ) ) μ 2 Ω 2 ( p ) 2 Φ 2 ( τ , p ) τ 2 + ε 2 Ω ( p ) [ Φ 2 2 ( τ , p ) λ ] Φ 2 ( τ , p ) τ + η 2 Φ 2 ( τ , p ) f 2 ( Φ 1 ( τ , p ) , Φ 2 ( τ , p ) , Φ 3 ( τ , p ) ) μ 3 Ω 2 ( p ) 2 Φ 3 ( τ , p ) τ 2 + ε 3 Ω ( p ) [ Φ 3 2 ( τ , p ) λ ] Φ 3 ( τ , p ) τ + η 3 Φ 3 ( τ , p ) f 3 ( Φ 1 ( τ , p ) , Φ 2 ( τ , p ) , Φ 3 ( τ , p ) ) )

不妨令 H ( t ) = 1 ,则可得零阶形变方程为

( 1 p ) L ( Φ 1 ( τ , p ) u 1 , 0 ( τ ) Φ 2 ( τ , p ) u 2 , 0 ( τ ) Φ 3 ( τ , p ) u 3 , 0 ( τ ) ) = p h N ( Φ 1 ( τ , p ) Φ 2 ( τ , p ) Φ 3 ( τ , p ) ) , (3-9)

上式满足初始条件

Φ 1 ( 0 , p ) = a ( p ) , Φ 1 ( τ , p ) τ | τ = 0 = b ( p ) , Φ 2 ( 0 , p ) = c ( p ) , Φ 2 ( τ , p ) τ | τ = 0 = d ( p ) , Φ 3 ( 0 , p ) = e ( p ) , Φ 3 ( τ , p ) τ | τ = 0 = 0. (3-10)

显然,我们很容易得到当 p = 0 时,方程(3-9)满足初始条件(3-10)的解为

Φ 1 ( τ , 0 ) = u 1 , 0 ( τ ) , Φ 2 ( τ , 0 ) = u 2 , 0 ( τ ) , Φ 3 ( τ , 0 ) = u 3 , 0 ( τ ) , Ω ( 0 ) = ω 0 . (3-11)

而当 p = 1 时,满足初始条件(3-10)的零阶形变方程(3-9)与满足初始条件(3-4)的系统(3-3)是等价的,因而有

Φ 1 ( τ , 1 ) = u 1 ( τ ) , Φ 2 ( τ , 1 ) = u 2 ( τ ) , Φ 3 ( τ , 1 ) = u 3 ( τ ) , Ω ( 1 ) = ω . (3-12)

由此可见当嵌入变量p从0变化到1时, Φ i ( τ , p ) 从初始猜测解 u i , 0 ( τ ) 变化到未知的精确解 u i ( τ ) ( i = 1 , 2 , 3 ) 。类似地,我们也可以得出 Ω ( p ) 从初始猜测频率 ω 0 往物理变化频率 ω 变化的结论。

我们把 Φ 1 ( τ , p ) , Φ 2 ( τ , p ) , Φ 3 ( τ , p ) Ω ( p ) 都展开成关于p的麦克劳林级数,可得到

Φ 1 ( τ , p ) = u 1 , 0 ( τ ) + m = 1 u 1 , m ( τ ) p m , Φ 2 ( τ , p ) = u 2 , 0 ( τ ) + m = 1 u 2 , m ( τ ) p m , Φ 3 ( τ , p ) = u 3 , 0 ( τ ) + m = 1 u 3 , m ( τ ) p m , Ω ( p ) = ω 0 + m = 1 ω m p m (3-13)

其中

u 1 , m ( τ ) = 1 m ! m Φ 1 ( τ , p ) p m | p = 0 , u 2 , m ( τ ) = 1 m ! m Φ 2 ( τ , p ) p m | p = 0 , u 3 , m ( τ ) = 1 m ! m Φ 3 ( τ , p ) p m | p = 0 , ω m = 1 m ! m Ω ( p ) p m | p = 0 . (3-14)

通过选取适当的h,我们可以使(3-13)式中的幂级数在 p = 1 时收敛。那么,结合(3-12)式与(3-13)式,我们可以得到方程的级数解为

u 1 ( τ ) = u 1 , 0 ( τ ) + m = 1 u 1 , m ( τ ) , u 2 ( τ ) = u 2 , 0 ( τ ) + m = 1 u 2 , m ( τ ) , u 3 ( τ ) = u 3 , 0 ( τ ) + m = 1 u 3 , m ( τ ) , ω = ω 0 + m = 1 ω m . (3-15)

u 1 , n = { u 1 , 0 ( τ ) , u 1 , 1 ( τ ) , , u 1 , n ( τ ) } , u 2 , n = { u 2 , 0 ( τ ) , u 2 , 1 ( τ ) , , u 2 , n ( τ ) } , u 3 , n = { u 3 , 0 ( τ ) , u 3 , 1 ( τ ) , , u 3 , n ( τ ) } , ω n = { ω 0 , ω 1 , , ω n } . (3-16)

同样地,我们对零阶形变方程(3-9)关于嵌入变量p求m次导数,再将求导得到的式子左右两边同除以 m ! ,并令 p = 0 ,则我们就可得到m阶形变方程

L ( u 1 , m ( τ ) χ m u 1 , m 1 ( τ ) u 2 , m ( τ ) χ m u 2 , m 1 ( τ ) u 3 , m ( τ ) χ m u 3 , m 1 ( τ ) ) = h ( R 1 , m ( u 1 , m 1 , ω m 1 ) R 2 , m ( u 2 , m 1 , ω m 1 ) R 3 , m ( u 3 , m 1 , ω m 1 ) ) , (3-17)

其中

( R 1 , m ( u 1 , m 1 , ω m 1 ) R 2 , m ( u 2 , m 1 , ω m 1 ) R 3 , m ( u 3 , m 1 , ω m 1 ) ) = 1 ( m 1 ) ! m 1 N ( Φ 1 ( τ , p ) Φ 2 ( τ , p ) Φ 3 ( τ , p ) ) p m 1 | p = 0 , (3-18)

此m阶形变方程满足以下初始条件

u 1 , m ( 0 ) = a m , u 1 , m ( 0 ) = b m , u 2 , m ( 0 ) = c m , u 2 , m ( 0 ) = d m , u 3 , m ( 0 ) = e m , u 3 , m ( 0 ) = 0 , ( m 1 ) (3-19)

由线性算子L的定义与解的表达原则可知,方程(3-17)的等式右边项里不能含有 cos τ sin τ ,也不能含有久期项 τ cos τ τ sin τ

我们把这些项都找出来,且令这些项的系数都为0,则可得

1 π 0 2 π R 1 , m ( u 1 , m 1 , ω m 1 ) cos τ d τ = 0 , 1 π 0 2 π R 1 , m ( u 1 , m 1 , ω m 1 ) sin τ d τ = 0 , 1 π 0 2 π R 2 , m ( u 2 , m 1 , ω m 1 ) cos τ d τ = 0 , 1 π 0 2 π R 2 , m ( u 2 , m 1 , ω m 1 ) sin τ d τ = 0 , 1 π 0 2 π R 3 , m ( u 3 , m 1 , ω m 1 ) cos τ d τ = 0 , 1 π 0 2 π R 3 , m ( u 3 , m 1 , ω m 1 ) sin τ d τ = 0 , (3-20)

其中

R 1 , m ( u 1 , m 1 , ω m 1 ) = n = 1 φ 1 ( m ) [ a m , n 1 cos n τ + b m , n 1 sin n τ ] , R 2 , m ( u 2 , m 1 , ω m 1 ) = n = 1 φ 2 ( m ) [ a m , n 2 cos n τ + b m , n 2 sin n τ ] , R 3 , m ( u 3 , m 1 , ω m 1 ) = n = 1 φ 2 ( m ) [ a m , n 3 cos n τ + b m , n 3 sin n τ ] . (3-21)

a m , n 1 = 1 π 0 2 π R 1 , m ( u 1 , m 1 , ω m 1 ) sin ( n τ ) d τ , b m , n 1 = 1 π 0 2 π R 1 , m ( u 1 , m 1 , ω m 1 ) sin ( n τ ) d τ , a m , n 2 = 1 π 0 2 π R 2 , m ( u 2 , m 1 , ω m 1 ) cos ( n τ ) d τ , b m , n 2 = 1 π 0 2 π R 2 , m ( u 2 , m 1 , ω m 1 ) sin ( n τ ) d τ , a m , n 3 = 1 π 0 2 π R 3 , m ( u 3 , m 1 , ω m 1 ) cos ( n τ ) d τ , b m , n 3 = 1 π 0 2 π R 3 , m ( u 3 , m 1 , ω m 1 ) sin ( n τ ) d τ , (3-22)

根据(3-17)和(3-20)式,我们可以求出 ω k , a k , b k , c k , d k , e k ( k = 0 , 1 , 2 , ) 的值。为使解的精度更高,我们需要对w进行一定程度的修正

ω = ω + v .

其中v为小参数。

4. 数值模拟与讨论

为了验证同伦分析方法应用于非线性多自由度系统的正确性和有效性,针对耦合van der Pol系统(3-2),我们将分四种情形对其数值模拟进行讨论,并对分别对这四种情形下的理论结果与数值结果进行比较。而耦合项 f i ( i = 1 , 2 , 3 ) 的取值可借鉴文献 [29] [30]。

4.1. 所有振子都同步运动

f 1 ( x 1 , x 2 , x 3 ) = x 1 x 2 f 2 ( x 1 , x 2 , x 3 ) = x 2 x 3 f 3 ( x 1 , x 2 , x 3 ) = x 3 x 1 (4-1)

我们取 μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 , 3 ) ,且令 a , b , c , d , e , ω 的初始猜测值分别为 a 0 = 2 , b 0 = 0 , c 0 = 2 , d 0 = 2 , e 0 = 2 , ω 0 = 1 ,我们需要选取合适的收敛控制参数h来确保级数解的收敛。为此,我们给出了图2所示的ω-h曲线来确定h的取值范围。细观图2可知,当 1.85 < h < 0.25 时,级数解 ω 是收敛的。为简便起见,我们取 h = 1 图3图4分别给出了八阶近似相图和时间历程曲线。

Figure 2.The eighth-order approximation of versus is obtained for μ i = 1 , ε i = 1 , η i = 1 ( i = 1 , 2 , 3 ) and λ = 1

图2. μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 , 3 ) 时的八阶ω-h近似曲线

Figure 3. Comparison of the phase curves of the eighth-order approximation with the numerical integration solution is given for μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 , 3 ) , a 0 = 2 , b 0 = 0 , c 0 = 2 , d 0 = 2 , e 0 = 2 and ω 0 = 1

图3. μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 , 3 ) , a 0 = 2 , b 0 = 0 , c 0 = 2 , d 0 = 2 , e 0 = 2 , ω 0 = 1 时的八阶近似相图与数值积分法的比较

Figure 4. Comparison of the time history responses of the eighth-order approximation with the numerical integration solution is given when μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 , 3 ) , a 0 = 2 , b 0 = 0 , c 0 = 2 , d 0 = 2 , e 0 = 2 and ω 0 = 1

图4. μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 , 3 ) , a 0 = 2 , b 0 = 0 , c 0 = 2 , d 0 = 2 , e 0 = 2 , ω 0 = 1 时的八阶时间历程曲线

利用相图(图3)与时间历程图(图4)可证明同伦分析方法求得的相图与数值分析方法所求得的相图结果比较吻合,且能明显看出 x 1 ( τ ) , x 2 ( τ ) , x 3 ( τ ) 三个振子同步运动。其中初始条件为:

x 1 ( 0 ) = 2.008600 , x 1 ( 0 ) = 0 , x 2 ( 0 ) = 2.008600 , x 2 ( 0 ) = 0 , x 3 ( 0 ) = 2.008600 , x 3 ( 0 ) = 0 ,

根据文献 [31] 并且可求得 x 1 ( τ ) , x 2 ( τ ) , x 3 ( τ ) , ω 的八阶近似表达式为

x 1 ( τ ) = 1.881842 cos ω ¯ τ + 0.720593 sin ω ¯ τ + 0.720593 cos 3 ω ¯ τ 1.881842 sin 3 ω ¯ τ + 0.158456 cos 5 ω ¯ τ 0.530028 sin 5 ω ¯ τ 0.176676 cos 7 ω ¯ τ 0.475368 sin 7 ω ¯ τ 0.020655 cos 9 ω ¯ τ 0.218080 sin 9 ω ¯ τ 0.043616 cos 11 ω ¯ τ + 0.103273 sin 11 ω ¯ τ 0.010428 cos 13 ω ¯ τ 0.000519 sin 13 ω ¯ τ 0.000074 cos 15 ω ¯ τ + 0.072994 sin 15 ω ¯ τ ,

x 2 ( τ ) = 1.881842 cos ω ¯ τ + 0.720593 sin ω ¯ τ + 0.720593 cos 3 ω ¯ τ 1.881842 sin 3 ω ¯ τ + 0.158456 cos 5 ω ¯ τ 0.530028 sin 5 ω ¯ τ 0.176676 cos 7 ω ¯ τ 0.475368 sin 7 ω ¯ τ 0.020655 cos 9 ω ¯ τ 0.218080 sin 9 ω ¯ τ 0.043616 cos 11 ω ¯ τ + 0.103273 sin 11 ω ¯ τ 0.010428 cos 13 ω ¯ τ 0.000519 sin 13 ω ¯ τ 0.000074 cos 15 ω ¯ τ + 0.072994 sin 15 ω ¯ τ ,

x 3 ( τ ) = 1.881842 cos ω ¯ τ + 0.720593 sin ω ¯ τ + 0.720593 cos 3 ω ¯ τ 1.881842 sin 3 ω ¯ τ + 0.158456 cos 5 ω ¯ τ 0.530028 sin 5 ω ¯ τ 0.176676 cos 7 ω ¯ τ 0.475368 sin 7 ω ¯ τ 0.020655 cos 9 ω ¯ τ 0.218080 sin 9 ω ¯ τ 0.043616 cos 11 ω ¯ τ + 0.103273 sin 11 ω ¯ τ 0.010428 cos 13 ω ¯ τ 0.000519 sin 13 ω ¯ τ 0.000074 cos 15 ω ¯ τ + 0.072994 sin 15 ω ¯ τ ,

ω = 0.942953.

4.2. 两个振子同步运动,而第三个振子以一无关的方式运动(除它与第二个振子有相同周期地振动外)

f 1 ( x 1 , x 2 , x 3 ) = x 1 x 2 f 2 ( x 1 , x 2 , x 3 ) = x 2 x 1 f 3 ( x 1 , x 2 , x 3 ) = x 1 x 2 10 x 3 (4-2)

在此情况中,我们取和情形4.1一样的参数值,即 μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 , 3 ) a , b , c , d , e , ω 的初始猜测值取为 a 0 = 2 , b 0 = 0 , c 0 = 2 , d 0 = 0 , e 0 = 2 , ω 0 = 1 ,我们可以从图5中确定出h的取值范围来确保级数解的收敛。因为h的取值范围为 1.05 < h < 0.65 ,为简便起见故仍然取 h = 1 ,一阶近似相图与时间历程图如图6图7所示。

Figure 5. The first-order approximation of ω versus h is obtained for μ i = 1 , ε i = 1 , η i = 1 ( i = 1 , 2 , 3 ) and λ = 1

图5. μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 , 3 ) 时一阶ω-h近似曲线

Figure 6. Comparison of the phase curves of the first-order approximation with the numerical integration solution is given for μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 , 3 ) , a 0 = 2 , b 0 = 0 , c 0 = 2 , d 0 = 2 , e 0 = 2 and ω 0 = 1

图6. μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 , 3 ) , a 0 = 2 , b 0 = 0 , c 0 = 2 , d 0 = 0 , e 0 = 2 , ω 0 = 1 时的一阶近似相图与数值积分法的比较

Figure 7. Comparison of the time history responses of the first-order approximation with the numerical integration solution is given when μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 , 3 ) , a 0 = 2 , b 0 = 0 , c 0 = 2 , d 0 = 2 , e 0 = 2 and ω 0 = 1

图7. μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 , 3 ) , a 0 = 2 , b 0 = 0 , c 0 = 2 , d 0 = 0 , e 0 = 2 , ω 0 = 1 时的一阶时间历程曲线

利用相图(图6)与时间历程图(图7)可证明同伦分析方法求得的相图与数值分析方法所求得的相图结果比较吻合。观察图7可发现,振子 x 1 ( τ ) x 2 ( τ ) 同步运动,而第三个振子 x 3 ( τ ) 则以一种无关于前两个振子的方式运动。其初始条件为 x 1 ( 0 ) = 4 , x 1 ( 0 ) = 0 , x 2 ( 0 ) = 4 , x 2 ( 0 ) = 0 , x 3 ( 0 ) = 4 , x 3 ( 0 ) = 0 ,且可得 x 1 ( τ ) , x 2 ( τ ) , x 3 ( τ ) , ω 近似表达式为:

x 1 ( τ ) = 1.844822 cos ω ¯ τ 0.118627 sin ω ¯ τ 0.273745 cos 3 ω ¯ τ 1.844822 sin 3 ω ¯ τ + 0.453207 cos 5 ω ¯ τ 0.120863 sin 5 ω ¯ τ 1.836718 cos 7 ω ¯ τ 0.135962 sin 7 ω ¯ τ 0.047096 cos 9 ω ¯ τ 0.125599 sin 9 ω ¯ τ 0.144917 cos 11 ω ¯ τ + 0.235481 sin 11 ω ¯ τ + 0.500366 cos 13 ω ¯ τ 0.003442 sin 13 ω ¯ τ 0.021553 cos 15 ω ¯ τ 0.035026 sin 15 ω ¯ τ ,

x 2 ( τ ) = 1.844822 cos ω ¯ τ 0.118627 sin ω ¯ τ 0.273745 cos 3 ω ¯ τ 1.844822 sin 3 ω ¯ τ + 0.453207 cos 5 ω ¯ τ 0.120863 sin 5 ω ¯ τ 1.836718 cos 7 ω ¯ τ 0.135962 sin 7 ω ¯ τ 0.047096 cos 9 ω ¯ τ 0.125599 sin 9 ω ¯ τ 0.144917 cos 11 ω ¯ τ + 0.235481 sin 11 ω ¯ τ + 0.500366 cos 13 ω ¯ τ 0.003442 sin 13 ω ¯ τ 0.021553 cos 15 ω ¯ τ 0.035026 sin 15 ω ¯ τ ,

x 3 ( τ ) = 1.124328 cos ω ¯ τ + 0.573061 sin ω ¯ τ + 0.134443 cos 3 ω ¯ τ 1.131251 sin 3 ω ¯ τ + 0.359854 cos 5 ω ¯ τ 0.670447 sin 5 ω ¯ τ 0.001711 cos 7 ω ¯ τ 0.223488 sin 7 ω ¯ τ + 0.271384 cos 9 ω ¯ τ + 0.348288 sin 9 ω ¯ τ + 0.054345 cos 11 ω ¯ τ 0.059881 sin 11 ω ¯ τ 0.171917 cos 13 ω ¯ τ 0.066372 sin 13 ω ¯ τ + 0.167520 cos 15 ω ¯ τ 0.089623 sin 15 ω ¯ τ ,

ω = 1.208175.

4.3. 环上相邻的振子之间彼此都相位差1/3周期的运动

f 1 ( x 1 , x 2 , x 3 ) = 3 / 9 ( x 1 + 2 x 2 ) 3 f 2 ( x 1 , x 2 , x 3 ) = 3 / 9 ( x 2 + 2 x 3 ) 3 f 3 ( x 1 , x 2 , x 3 ) = 3 / 9 ( x 3 + 2 x 1 ) 3 (4-3)

在此情况中,环上相邻的振子之间彼此都相位差 周期的运动,因而对系统(3-1)中的各系数取值与前两种情况有所不同。这里我们取 μ i = 1 , ε i = 1 , η i = 4 , λ = 1 ( i = 1 , 2 , 3 ) 。同样地,初始猜测值也被赋予一些新的数值, a 0 = 10 / 5 , b 0 = 30 / 5 , c 0 = 10 / 5 , d 0 = 30 / 5 , e 0 = 2 10 / 5 , ω 0 = 2 。为确保级数解的收敛,我们取 h = 1 ,从而得到相应地相图与时间历程图(见图8图9)。

Figure 8. Comparison of the phase curves of the first-order approximation with the numerical integration solution is given for μ i = 1 , ε i = 1 , η i = 4 , λ = 1 ( i = 1 , 2 , 3 ) , a 0 = 10 / 5 , b 0 = 30 / 5 , c 0 = 10 / 5 , d 0 = 30 / 5 , e 0 = 2 10 / 5 and ω 0 = 2

图8. μ i = 1 , ε i = 1 , η i = 4 , λ = 1 ( i = 1 , 2 , 3 ) , a 0 = 10 / 5 , b 0 = 30 / 5 , c 0 = 10 / 5 , d 0 = 30 / 5 , e 0 = 2 10 / 5 , ω 0 = 2 时的一阶近似相图与数值积分法的比较

Figure 9. Comparison of the time history responses of the first-order approximation with the numerical integration solution is given for μ i = 1 , ε i = 1 , η i = 4 , λ = 1 ( i = 1 , 2 , 3 ) , a 0 = 10 / 5 , b 0 = 30 / 5 , c 0 = 10 / 5 , d 0 = 30 / 5 , e 0 = 2 10 / 5 and ω 0 = 2

图9. μ i = 1 , ε i = 1 , η i = 4 , λ = 1 ( i = 1 , 2 , 3 ) , a 0 = 10 / 5 , b 0 = 30 / 5 , c 0 = 10 / 5 , d 0 = 30 / 5 , e 0 = 2 10 / 5 , ω 0 = 2 时的一阶时间历程曲线

Figure 10. Comparison of time history diagrams of x 1 , x 2 and x 3 in the same condition

图10. x 1 , x 2 , x 3 在同一条件下时间历程图的比较

图8可知,利用同伦分析方法求得的近似解与数值积分法求得的解基本吻合,而由图10也可看出三个振子正相差1/3周期运动。我们同样可求得数值积分法的初始条件为

x 1 ( 0 ) = 0.632455 , x 1 ( 0 ) = 1.095445 , x 2 ( 0 ) = 0.632455 , x 2 ( 0 ) = 1.095445 , x 3 ( 0 ) = 1.264911 , x 3 ( 0 ) = 0

我们也就能得到 x 1 ( τ ) , x 2 ( τ ) , x 3 ( τ ) ω 的近似解表达式:

x 1 ( τ ) = 0.649572 cos τ + 1.089516 sin τ + 0.017116 cos 3 τ + 0.001976 sin 3 τ ,

x 2 ( τ ) = 0.615339 cos τ 1.101374 sin τ 0.017116 cos 3 τ + 0.001976 sin 3 τ ,

x 3 ( τ ) = 1.264911 cos τ + 0.094868 sin τ 0.031623 sin 3 τ ,

ω = 2.

4.4. 两个振子相位差1/2周期,而第三个振子2倍于它们的频率振动

f 1 ( x 1 , x 2 , x 3 ) = ( 1 x 1 2 ) x 2 f 2 ( x 1 , x 2 , x 3 ) = ( 1 x 2 2 ) x 1 f 3 ( x 1 , x 2 , x 3 ) = x 1 3 x 2 3 (4-4)

此时,两个振子相位差1/2周期,而第三个振子 倍于它们的频率振动。为计算方便,我们取 μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 ) , μ 3 = 1 , ε 3 = 1 / 16 , η 3 = 4

而初始猜测值求得为 a 0 = 1 , b 0 = 0 , c 0 = 1 , d 0 = 0 , e 0 = 1 , ω 0 = 1 。取 h = 1 得近似解的相图和时间历程图如图11图12所示。

Figure 11. Comparison of the phase curves of the first-order approximation with the numerical integration solution is given for μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 ) , μ 3 = 1 , ε 3 = 1 / 16 , η 3 = 4 , a 0 = 1 , b 0 = 0 , c 0 = 1 , d 0 = 0 , e 0 = 1 , ω 0 = 1

图11. μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 ) , μ 3 = 1 , ε 3 = 1 / 16 , η 3 = 4 , a 0 = 1 , b 0 = 0 , c 0 = 1 , d 0 = 0 , e 0 = 1 , ω 0 = 1 时的一阶近似相图

Figure 12. Comparison of the time history responses of the first-order approximation with the numerical integration solution is given for μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 ) , μ 3 = 1 , ε 3 = 1 / 16 , η 3 = 4 , a 0 = 1 , b 0 = 0 , c 0 = 1 , d 0 = 0 , e 0 = 1 and ω 0 = 1

图12. μ i = 1 , ε i = 1 , η i = 1 , λ = 1 ( i = 1 , 2 ) , μ 3 = 1 , ε 3 = 1 / 16 , η 3 = 4 , a 0 = 1 , b 0 = 0 , c 0 = 1 , d 0 = 0 , e 0 = 1 , ω 0 = 1 时的一阶近似时间历程曲线

Figure 13. Comparison of time history diagrams of x 1 , x 2 and x 3 in the same condition

图13. x 1 , x 2 , x 3 在同一条件下时间历程图的比较

通过观察与计算可知,利用同伦法得到的相图与时间历程曲线是合适的,再次证明了同伦分析方法于耦合多自由度振子环系统上的适用性。其数值积分法满足的初始条件为

x 1 ( 0 ) = 1 , x 1 ( 0 ) = 0 , x 2 ( 0 ) = 1 , x 2 ( 0 ) = 0 , x 3 ( 0 ) = 1 , x 3 ( 0 ) = 0

同样地,得到 x 1 ( τ ) , x 2 ( τ ) , x 3 ( τ ) 和频率 ω 的一阶近似解为:

x 1 ( τ ) = cos τ ,

x 2 ( τ ) = cos τ ,

x 3 ( τ ) = cos 2 τ ,

ω = 1.

显然,此近似解就是满足条件(4-4)的系统(3-1)的精确解,也符合第四种情况的振子运动(见图13)。

5. 结论

本文利用该方法对多自由度耦合van der Pol振子环进行了精确的近似分析。由之前的数值模拟可以知道,利用同伦分析方法求解三自由度耦合van der Pol方程得到的周期解比较精确。该方法不仅可以计算任意阶次的频率和位移,而且计算量小。与经典摄动技术相比,该方法不需要依靠小参数并可以很好地保证级数解的收敛性以及解的精度,在解决非线性问题方面具有非常好的效果。在后续研究中,我们可以考虑四种不同情况中的振子如何才能互相转换,在何种外激励下使系统形成四稳态局面。

致谢

感谢国家自然科学基金项目(11572288)、浙江省自然科学基金项目(LY20A020003)和浙江师范大学数计学院重点项目(ZD01201713)资助。

NOTES

*通讯作者。

文章引用:
胡佩倩, 钱有华. 一类环上的三自由度耦合van der Pol方程的同伦分析方法[J]. 动力系统与控制, 2020, 9(1): 22-39. https://doi.org/10.12677/DSC.2020.91003

参考文献

[1] Chen, H., Huang, J.L. and Sze, K.Y. (2007) Multidimensional Lindstedt-Poincare Method for Nonlinear Vibration of Axially Moving Beams. Journal of Sound and Vibration, 306, 1-11.
https://doi.org/10.1016/j.jsv.2007.05.038
[2] Liao, S.J. (2009) Notes on the Homotopy Analysis Method: Some Definitions and Theorems. Communication in Nonlinear Sciences and Numerical Simulation, 14, 983-997.
https://doi.org/10.1016/j.cnsns.2008.04.013
[3] Marinca, V. and Herisanu, N. (2010) Determination of Periodic Solutions for the Motion of a Particle on a Rotating Parabola by Means of the Optimal Homotopy Asymptotic Method. Journal of Sound and Vibration, 329, 1450-1459.
https://doi.org/10.1016/j.jsv.2009.11.005
[4] Liao, S.J. (2010) An Optimal Homotopy-Analysis Approach for Strongly Nonlinear Differential Equations. Communication in Nonlinear Sciences and Numerical Simulation, 15, 2003-2016.
https://doi.org/10.1016/j.cnsns.2009.09.002
[5] Niu, Z. and Wang, C. (2010) A One-Step Optimal Homotopy Analysis Method for Nonlinear Differential Equations. Communication in Nonlinear Sciences and Numerical Simulation, 15, 2026-2036.
https://doi.org/10.1016/j.cnsns.2009.08.014
[6] Chen, S.H., Chen, Y.Y. and Sze, K.Y. (2009) A Hyperbolic Perturbation Method for Determining Homoclinic Solution of Certain Strongly Nonlinear Autonomous Oscillators. Journal of Sound and Vibration, 322, 381-392.
https://doi.org/10.1016/j.jsv.2008.11.015
[7] 廖世俊. 超越摄动: 同伦分析方法基本思想及其应用[J]. 力学进展, 2008, 38(1): 1-34.
[8] 廖世俊. 同伦分析方法: 一种不依赖于小参数的非线性分析方法[J]. 上海力学, 1997, 18(3): 196-201.
[9] Liao, S.J. and Chwang, A.T. (1998) Application of Homotopy Analysis Method in Nonlinear Oscillation. ASME Journal of Application Mechanics, 65, 914-922.
https://doi.org/10.1115/1.2791935
[10] Rana, J. and Liao, S.J. (2019) A General Analytical Approach to Study Solute Dispersion in Non-Newtonian Fluidflow. European Journal of Mechanics B: Fluids, 77, 183-200.
https://doi.org/10.1016/j.euromechflu.2019.04.013
[11] Qian, Y.H., Zhang, W., Lin, B.W. and Lai, S.K. (2011) Analytical Approximate Periodic Solutions of Multi-Degree-of-Freedom Strongly Nonlinear Coupled van der Pol Oscillators by Homotopy Analysis Method. Acta Mechanica, 219, 1-14.
https://doi.org/10.1007/s00707-010-0433-3
[12] 崔继峰. 受迫Van der Pol-Duffing振子的倍周期与拟周期之解析解研究[C]//中国力学学会动力学与控制专业委员会. 第十届动力学与控制学术会议摘要集. 中国力学学会动力学与控制专业委员会: 中国力学学会, 2016: 276.
[13] 许天亮, 樊晓敏, 张跃进. 非线性分数阶微分方程的同伦分析解法[J]. 湘潭大学自然科学学报, 2016, 38(4): 6-9.
[14] 郭江明. 时滞耦合van der Pol系统的双Hopf分岔与周期解的同伦方法研究[D]: [硕士学位论文]. 金华: 浙江师范大学, 2017.
[15] 杨孟, 林炎海. 圆盘上非牛顿流体Marangoni边界层同伦分析解研究[C]//第十届全国流体力学学术会议论文摘要集. 华侨大学数学科学学院计算机科学福建省高校重点实验室, 2018: 121.
[16] Rompala, K., Rand, R and Howland, H. (2007) Dynamics of Three Coupled van der Pol Oscillators with Application to Circadian Rhythms. Communications in Nonlinear Science and Numerical Simulation, 12, 794-803.
https://doi.org/10.1016/j.cnsns.2005.08.002
[17] Tang, J.S. and Han, F. (2009) Amplitude Control of a Limit Cycle in a Coupled van der Pol System. Nonlinear Analysis: Theory, Methods & Applications, 7, 2491-2496.
https://doi.org/10.1016/j.na.2009.01.130
[18] 钱有华, 张伟, 姚志刚. 参数激励非自治薄板系统周期解的同伦分析方法[C]//第十二届全国非线性振动暨第九届全国非线性动力学和运动稳定性学术会议论文集. 2009.
[19] 吴俊. 两个耦合Van der Pol方程拟周期解的存在性[D]: [硕士学位论文]. 长沙: 湖南师范大学, 2016.
[20] 石兆羽, 杨绍普, 赵志宏. 基于耦合van der Pol-Duffing系统的微弱信号检测研究[J]. 中国测试, 2018, 44(8): 107-112.
[21] 施添添, 茅晓晨. 时滞耦合van der Pol-Duffing振子环的动力学分析[J]. 动力学与控制学报, 2019, 17(3): 264-269.
[22] Vinod, V., Bipin, B., Narayanan, M.D. and Mihir, S. (2015) Effect of Oscillator and Initial Condition Differences in the Dynamics of a Ring of Dissipative Coupled van der Pol Oscillators. Journal of Mechanical Science and Technology, 29, 1931-1939.
https://doi.org/10.1007/s12206-015-0103-4
[23] Barron, M.A. (2016) Stability of a Ring of Coupled van der Pol Oscillators with Non-Uniform Distribution of the Coupling Parameter. Journal of Applied Research and Technology, 14, 62-66.
https://doi.org/10.1016/j.jart.2016.01.002
[24] Fuad, R., Lu, B.L. and Surendra, S. (2007) Synchronized and Unsynchronized Chaos in a Modulated Bidirectional Ring Laser. Physical Review A: Atomic, Molecular and Optical Physics, 75, Article ID: 043811.
https://doi.org/10.1103/PhysRevA.75.043811
[25] Fernando, A., Ana, P.S.D. and Carla, M.A.P. (2010) Quasi-Periodic States in Coupled Rings of Cells. Communications in Nonlinear Science and Numerical Simulation, 15, 1048-1062.
https://doi.org/10.1016/j.cnsns.2009.05.035
[26] Faralli, S., Gambini, F., et al. (2016) Bidirectional Transmission in an Optical Network on Chip with Bus and Ring Topologies. IEEE Photonics Journal, 8, Article No. 0600407.
https://doi.org/10.1109/JPHOT.2016.2526607
[27] Xu, X., Yu, D.Y. and Wang, Z.H. (2019) Inter-Layer Synchronization of Periodic Solutions in Two Coupled Rings with Time Delay. Physica D: Nonlinear Phenomena, 396, 1-11.
https://doi.org/10.1016/j.physd.2019.02.010
[28] 张锁春. 可激励系统分析的数学理论[M]. 北京: 科学出版社, 2010.
[29] 段春梅. 一类耦合van der Pol系统的改进同伦分析方法及应用[D]: [硕士学位论文]. 金华: 浙江师范大学, 2012.
[30] Nohara, B.T. and Arimoto, A. (2009) Nonexistence Theorem Except the Out-of-Phase and in-Phase Solutions in the Coupled van der Pol Equation System. Ukrainian Mathematical Journal, 61, 1311-1337.
https://doi.org/10.1007/s11253-010-0278-x
[31] Qian, Y.H. and Chen, S.M. (2010) Accurate Approximate Analytical Solutions for Multi-Degree-of-Freedom Coupled van der Pol-Duffing Oscillators by Homotopy Analysis Method. Communications in Nonlinear Science and Numerical Simulation, 15, 3113-3130.
https://doi.org/10.1016/j.cnsns.2009.11.027