具非均匀剪切项的耦合振子系统的同步分析

李娜 刘辉昭

摘要 讨论了具非均匀剪切项的全局耦合相位振子系统由非同步状态向同步状态转变的动力学问题，其中剪切项的强度服从洛伦兹分布。利用改进的Kuramoto模型，通过运用OA约化方法，得到了相位振子系统的约化方程，进而利用微分方程的稳定性理论，得到了系统标准差阈值，这里得到标准差阈值方法与现有文献中通过构造拟设的方法所求得标准差阈值方法更为一般，此外还通过作出约化方程的波形图、相图以及序参量随时间变化的图像来进一步验证所得结果。当剪切强度标准差超过该阈值时，系统在任意耦合强度下不可能同步，当剪切强度标准差小于该阈值时，系统在耦合强度超过某一定值时会产生同步。這与在无剪切项作用下，系统可以通过调整耦合强度最终达到同步状态的结果不同，说明了剪切项对耦合相位振子系统具有重要的影响。

关 键 词 耦合振子系统；Kuramoto模型；剪切项；标准差阈值；同步

中图分类号 O29 文献标志码 A

同步现象涉及物理、生物、化学等自然科学、工程学甚至社会科学等许多领域[1-6]。该类问题最早是由惠更斯于1673年发现2个相邻的钟摆出现同步摆动而引发学者广泛关注。1955年，Wiener[7]开始在数学上来研究含大量耦合振子的系统同步现象。1967年Winfree[8]利用平均场理论给出了简化的耦合振子摸型。1975年，Kuramoto[9]在Winfree的研究基础上提出了如式（1）的全局耦合相位振子平均场模型：

[θi（t）=ωi+KNj=1Nsin（θj（t）-θi（t）），i=1，2，...，N，] （1）

式中：[θi]为第[i]个振子的相位；[ωi]为第[i]个振子的满足一定分布的固有频率；[N]为系统总振子数目；[K]为耦合强度。

为了刻画系统中振子相位的有序程度，通常引入如式（2）的复序参量：

[r（t）=R（t）eiψ（t）=N-1k=1Neiθk（t）]， （2）

来描述系统的有序程度，式中：[ψ（t）]为系统中所有振子的平均相位；[R（t）]表示系统中所有振子的平均振幅。文献[10]指出，当系统达到完全同步时，即[θi（t）=θj（t），∀i，j=1，2，...，N]，此时[r（t）=1]；而当振子之间相位没有固定关系时（非同步），在[N→∞]时，[r（t）=0]，因此序参量可以描述振子系统有序程度。

Kuramoto[11]利用平均场理论得到[r（t）]的自洽方程，进而得到如下的结论：当耦合强度从零增加至某一定值[K0]时，序参量会经历从零到非零转变，相应的系统随着耦合强度的增加由非同步状态转变为同步状态。

近年来关于耦合振子系统的同步研究持续增多，Tanaka等考虑了系统中有阻尼作用的情况[12-13]，在原模型中添加了惯性项，发现在具有阻尼作用下系统出现不连续的相变和滞后现象，与实验结果相吻合。文献[14]研究了振子系统在单向最近邻耦合项作用下的情况，利用改进简化的Kuramoto相位模型，发现振子数在小于临界值时系统只出现同步状态或者非同步状态，不会出现部分同步状态，而当振子数目大于该临界值时，系统呈现规律的同步状态与非同步状态共存的现象，并且对系统中出现的部分同步状态规律及其稳定性进行了理论分析，得到了部分同步状态的渐近稳定解。文献[15]讨论了耦合振子系统在具全局耦合作用的网络上，分别在有噪声和无噪声的情况下的同步现象，并且得到了相应的振子系统平均场的振幅与频率的显式解。

剪切项是振子振动出现复杂行为的非常关键的非线性成分[9]，亦是振子系统重要的组成成分，剪切项对耦合振子系统的作用却是不可忽视的，但在耦合振子系统中关于剪切项的研究是稀缺的，所以对具剪切项的耦合振子系统的研究是非常有必要的。在剪切项固定（为常数）时，文献[17]通过对实际的耦合纳米机械振子组成的耦合振子系统的研究，得到系统达到同步（锁相）状态的边界解，并得到系统从非同步态向同步状态转化的结果。

文献[18]通过OA流形理论，并利用构造拟设的方法研究了耦合相位振子在非均匀剪切项的作用下振子族之间的同步动力学，得到了临界耦合强度与剪切强度的关系式，并说明了剪切项的标准差可能阻碍系统同步的发生。

本文主要研究振子之间的非均匀剪切项在服从指定分布时，对耦合振子系统的动力学行为的影响，利用OA流形约化方法以及微分方程定性理论，发现了在非均匀剪切项作用下，剪切项的标准差对耦合振子系统的同步行为有重要的影响，在超过剪切项标准差的阈值时，会阻止系统中振子同步行为的发生，耦合振子系统处于稳定的非同步状态。

本文中根据经典的Kuramoto相位模型考虑具非均匀剪切项的全局耦合相位振子系统，系统描述如下[9]，其中振子i所受系统内其他振子的非均匀剪切项表示为[qi]：

[θi（t）=ωi+Kqi+KNj=1N[sin（θj（t）-θi（t））-qicos（θj（t）-θi（t））]，i=1，2，...，N]。 （3）

1 理论分析

在系统式（3）中引入复序参量[r（t）]可得

[θi（t）=ωi+Kqi+KR（t）[sin（ψ（t）-θi（t））-qicos（ψ（t）-θi（t））]]。 （4）

考虑系统在[N→∞]的热力学极限情况下，[N]个耦合振子的运动由如下的方程描述[19]：

[θ=v（θ，t）]。 （5）

本文中所考虑的具剪切项的耦合振子系统式（3）中的速度项[v（θ，t）]为

[v（θ，t）=ω+Kq+K2i（r（t）e-iθ（1-iq）-r*（t）eiθ（1+iq））]。 （6）

引入相位的概率密度函数[fθ，ω，q，t][20]，则[fθ，ω，q，tdθdωdq]表示振子相位处于[θ]到[θ+dθ]，固有频率处于[ω]到[ω+dω]，剪切强度处于[q]到[q+dq]之间的振子比率。

[r（t）]是序参量在[t]时刻的值，可以用来刻画耦合系统在[t]时刻的状态，[r*（t）]表示[r（t）]的共轭

[r（t）=-∞∞-∞∞02πf（θ，ω，q，t）eiθdθdωdq]， （7）

由不存在扩散项的Feynman-Kac公式，得到系统式（5）等价于不带扩散项的Fokker-Planck方程，即连续方程

[∂∂tf+∂∂θ（vf）=0]。 （8）

由上述讨论可得密度函数[fθ，ω，q，t]服从连续方程

[∂∂tf+∂∂θ（{ω+Kq+K2i[r（t）e-iθ（1-iq）-c.c.]}f）=0]， （9）

式中，c.c.代表前一项[r（t）e-iθ（1-iq）]的复共轭量。

由于密度函数[fθ，ω，q，t]是实的且其中[θ]以[2π]为周期的，可将[fθ，ω，q，t]以傅里叶展开式展开，得到

[fθ，ω，q，t=p（ω，q）2πl=-∞∞fl（ω，q，t）eilθ]， （10）

式中：[fl=f*-l]；[f0=1]；且[p（ω，q）]是固有频率[ω]和剪切强度[q]的联合概率密度函数。密度函数[fθ，ω，q，t]的傅里叶展开式的第1项决定了序参量因而十分重要。

[r*（t）=-∞∞-∞∞p（ω，q）f1（ω，q，t）dωdq]。 （11）

将傅里叶展开式（10）代入到连续性方程（9）中可得到

[∂∂tf+il（ω+Kq）fl-Kl2[r*（t）（1+iq）fl-1-r（t）（1-iq）fl+1]=0]， （12）

由Ott和Antonsen发现下列拟设是Kuramoto模型带有固有频率分布系统式（3）的一个特殊解[21]

[fl（ω，q，t）=a（ω，q，t）l]。 （13）

由OA流形约化理论[21]，在此处借助于上面得到的积分微分方程（7）、（12），将[a（ω，q，t）]代入到方程（7）、（12）中得到

[∂a（ω，q，t）∂t+i（ω+kq）a（ω，q，t）+k2[r（t）（1-iq）a（ω，q，t）2-r*（t）（1+iq）]=0]。 （14）

在上述讨论中假设固有频率[ω]与剪切强度[q]是彼此相互独立的随机变量，即有[p（ω，q，t）=g（ω）h（q）]，在本文中取固有频率[ω]与剪切强度[q]分别具有如下的服从洛伦兹分布的概率密度函数：

[g（ω）=δπ（ω-ω0）2+δ2，h（q）=γπ（q-q0）2+γ2，] （15）

式中：[ω0]为固有频率的均值；[δ]为固有频率的标准差；[q0]为剪切强度的均值；[γ]为剪切强度标准差。因而[δ]以及[γ]均不小于0，即[δ≥0]和[γ≥0]。

在系统中的振子的固有频率[ω]以及剪切项的强度[q]分别服从的洛伦兹分布式（15）可以化为下列形式：

[g（ω）=12πi（1ω-ω0-iδ-1ω-ω0+iδ） ，h（q）=12πi（1q-q0-iγ-1q-q0+iγ） 。] （16）

再者，[a（ω，q，t）]在选取适当的[ω]复平面上是解析延拓的[22]，那么方程（11）可以通过在复[ω]平面的下半平面上作一无穷大的围道，利用留数定理可计算出关于固有频率[ω]的极点[ωp=ω0-iδ]的积分值

[r（t）=-∞∞a*（ωp，q，t）h（q）dq]。 （17）

同样的选取适当的[q]复平面，使[a（ωp，q，t）]在[q]复平面上是的解析延拓的，在[q]复平面的下半平面上利用留数定理计算出关于剪切项强度的极点[qp=q-iγ]的积分

[r（t）=a*（ωp，qp，t）]。 （18）

令式（18）中[ωp=ω0-iδ] ，[qp=q0-iγ]，并将得到结果代入式（14）中，可以得到具剪切项的耦合振子系统的约化方程

[r*（t）=-i（ω0-iδ+K（q0-iγ））r*（t）+K2{r*（t）[1+i（q0-iγ）]-r（t）[1-i（q0-iγ）]r*（t）2}]。 （19）

在接下来的讨论中，不妨考虑振子之间的耦合强度[K>0]。为了使对复序参量[r（t）]讨论更简便直观，令[r*（t）=x（t）+iy（t）]代入式（19）中可得

[x（t）+iy（t）=-（ω0-iδ+K（q0-iγ））（ix（t）-y（t））+K2{x（t）-（q0-iγ）y（t）+i[（q0-iγ）x（t）+y（t）]-]

[[x（t）-（q0-iγ）y（t）]（x2（t）-y2（t））-i[x（t）-（q0-iγ）y（t）]2x（t）y（t）+]

[i[（q0-iγ）x（t）+y（t）]（x2（t）-y2（t））-[（q0-iγ）x（t）+y（t）]2x（t）y（t）}。] （20）

分离式（20）的实虚部得到

[x（t）=-（δ+K2（γ-1））x（t）+（ω0+K2q0）y（t）+K2（γ-1）x3（t）-K2q0y3（t）-K2（1-γ）x（t）y2（t）-K2q0x2（t）y（t） ，y（t）=-（ω0+K2q0）x（t）-（δ+K2（γ-1））y（t）+K2q0x3（t）+K2（γ-1）y3（t）+K2q0x（t）y2（t）+K2（γ-1）x2（t）y（t） 。] （21）

微分方程組（21）就是以复序参量为变量式（19）的等价方程，亦是具剪切项的耦合振子系统式（3）的约化系统。通过研究等价方程组（21）的动力学，进而研究复序参量的演化，最终可以得到有关具剪切项的耦合振子系统式（3）的有关同步状态和非同步状态的动力学。

显然，（0， 0）是微分方程组（21）的平衡点。

再者令

[Φ（x，y）=K2（γ-1）x3（t）-K2q0y3（t）-K2（1-γ）x（t）y2（t）-K2q0x2（t）y（t） ，Ψ（x，y）=K2q0x3（t）+K2（γ-1）y3（t）+K2q0x（t）y2（t）-K2（1-γ）x2（t）y（t） 。] （22）

在（0， 0）邻域内，[Φ（x，y），Ψ（x，y）]对[x，y]连续可微，且[Φ（0，0）=Ψ（0，0）=0]，则方程（21）在平衡点（0， 0）处可以得到如下线性近似系统：

[x（t）=-（δ+K2（γ-1））x（t）+（ω0+K2q0）y（t） ，y（t）=-（ω0+K2q0）x（t）-（δ+K2（γ-1））y（t） 。] （23）

线性近似系统式（23）在（0， 0）处行列式为

[I0，0=-（δ+K2（γ-1））-（ω0+K2q0）ω0+K2q0-（δ+K2（γ-1））]。 （24）

可得到特征方程为

[（λ+δ+K2（γ-1））2+（ω0+K2q0）2=0]。 （25）

从而系统式（23）在（0， 0）点对应的特征值为

[λ1，2=-δ+K2（1-γ）±i（ω0+K2q0）]。

当[Re（λ1，2）<0]时，即[γ>1-2δK]时，[（0，0）]点是线性近似系统式（23）的稳定的平衡点，由于[（0，0）]是原系统式（21）的双曲平衡点[23]，[Φ（x，y），Ψ（x，y）]在原点邻域内对[x，y]连续可微，且[Φ（x，y），Ψ（x，y）=o（r），][r→0（r=x2+y2）]，因此[（0，0）]亦是原系统式（21）的稳定的平衡点。

当[Re（λ1，2）=0]，即[γ=1-2δK]时，对原系统式（21）进行极坐标变换化为

[r=K（γ-1）2r3 ，θ=-ω0+Kq02（r2-1） 。] （26）

对式（26）积分后可得到解为

[r=±1-K（γ-1）t-C1 ，θ=-（ω0+Kq02）t-lntC1K（γ-1）+C2 ，] （27）

式（27）中[C1，C2]为任意常数。由此可见，当[t→+∞]时，此时剪切强度标准差需满足[γ<1]以保证[r]有意义，[r→0，θ→+∞]，因此[（0，0）]是系统式（21）稳定的焦点。

当[Re（λ1，2）>0]，即[γ>1-2δK]时， 显然[（0，0）]在线性近似系统式（23）中是不稳定的，从而在系统式（21）中[（0，0）]是不稳定的。

通过以上分析，以及微分方程的定性理论[23]，对于系统式（21）的平衡点[（0，0）]有如下的结论：当[γ≤1-2δK]时，[（0，0）]是系统式（21）不稳定的平衡点；而当 [γ>1-2δK]时，[（0，0）]是系统式（21）稳定的平衡点。

最后，由前面复序参量的定义式（2）以及[r*（t）=x（t）+iy（t）]可知，在[（0，0）]点处，此时[x（t）=0，y（t）=0]，相应的复序参量的模[r（t）=0]，对应的耦合振子系统处于非同步状态。

根据以上分析可以得到，当[γ≤1-2δK]时，具剪切项的耦合振子系统处于不稳定的非同步状态，耦合振子系统可能会发生向同步状态的转化；当[γ>1-2δK] 时，耦合振子系统处于稳定的非同步状态，此时具剪切项的耦合振子系统在剪切项的作用下振子不能达到同步状态。故而具剪切项的耦合振子系统的剪切项的强度的标准差的临界值为[γc=1-2δK]。

下面，通过数值模拟来进一步验证本文中得到的结论。

2 数值模拟

对系统式（21）进行分析，选取参数值：[K=2，q0=0.5，ω0=3，δ=0.1]。利用上文理论分析可知剪切强度的标准差临界值[γc=1-2δK=1-0.1=0.9]。下面通过MATLAB作出系统式（21）的波形图与相图，来说明剪切项强度的标准差[γ]对耦合振子系统同步状态的影响，得到与上文分析的结果是相一致。

通过选取不同的剪切强度的标准差[γ]来对比剪切强度标准差对耦合振子系统的同步状态与非同步状态的稳定性的影响。分别取[γ=0]，[γ=0.85]，[γ=0.90]，[γ=0.95]得到结果如图1～图6。

根据前面参数选择可知晓[γc=1-2δK=0.9]，分别取[γ]的值为0，0.85，0.90，0.95，作出系统式（21）的波形图、相图以及复序参量的模值随时间变化图。从图1可以观察到，当[γ=0]时，此时系统中的剪切项的强度固定（[q]为常数），系统式（21）平衡点（0， 0）是不稳定且有周期解产生，相对应的是系统式（19）处于部分同步状态，再者由图2序参量的的变化趋势是在较短的时间内从0向1趋近，可知耦合振子系统从非同步态迅速地转化为稳定的同步态；随着[γ]从0增大到0.85，此时[γ<γc]，由图3可观察到系统式（21）的平衡点是不稳定的，系统中的周期解个数逐渐增多，相应的系统式（19）处于稳定的部分同步状态，由图4可知序参量的模值[r（t）]也是从零逐渐趋于某一常数（大于零），与图2相比，在相同的耦合强度的作用下，[γ=0]时序参量[r（t）]在8秒内达到稳定值0.97，而在[γ=0.85]时，序参量[r（t）]在约在100秒达到定值0.56，说明随着耦合振子系统内剪切强度的标准差的增大，耦合振子系统能够达到稳定的同步状态或部分同步状态所需时间逐渐增加，系统中振子能达到的同步的比例是逐渐减小的。当[γ]增大至0.90时，此时[γ=γc]，由图5可知系统式（21）的平衡点是不稳定的，系统中周期解的个数逐渐的减少到1个，由图6序参量的值[r（t）]接近于0，从0.012逐渐减小到0，说明在此剪切强度标准差作用下，耦合振子系统由不稳定的非同步状态向稳定的非同步态转变；当[γ]达到0.95时，即此时[γ>γc]，图7可观察到系统式（21）的平衡点是渐近稳定的，图8中序参量[r（t）]接近于1，且在一段時间内能较快速地趋于0，即对应的振子系统较快的达到稳定的非同步状态。由以上分析，可以说明在耦合强度一定时，剪切强度的标准差[γ]从零增大到[1-2δK]，耦合振子系统是处于同步或部分同步的状态，但超过该值时，耦合振子系统处于非同步状态，以上分析结果与前文中结论相一致。

本節中对具服从指定分布函数剪切项的耦合振子系统应用Matlab做出数值模拟结果，选取适当的参数值分别画出了系统的波形图与相图，研究剪切项标准差对耦合振子系统同步状态的影响，发现能较好的说明之前得到的关于具剪切项耦合阵子系统的剪切项强度标准差对系统同步有极大的影响，验证了存在剪切项的临界标准差，当耦合振子系统在超过这一临界值时，系统就不能达到同步状态。并且也作出复序参量在一定的剪切项标准差下随时间变化趋势的图像，进一步说明了在剪切强度的标准差超过临界值时，耦合振子系统在非均匀剪切项作用下不能达到同步状态。

3 结论

本文利用改进的Kuramoto模型研究了具非均匀剪切项对全局耦合相位振子系统同步行为的影响以及作用，通过引入序参量来描述系统的同步（有序）程度，得到关于序参量的微分方程组。通过对该方程组的稳定性分析，得到了系统稳定性与剪切强度标准差的关系，从而揭示了具剪切项Kuramoto系统的剪切强度对系统同步的影响。利用数值模拟，进一步证明了剪切强度的标准差在较大时会阻止耦合振子系统的同步，使耦合振子系统在充分大的耦合强度作用下保持非同步状态。本文方法更具一般性.本文研究的是具剪切项的全局耦合的Kuramoto模型，而在耦合结构不同时，例如局部耦合结构下，具剪切项的振子系统的同步动力学行为以及剪切项对耦合振子系统同步的影响等将在今后工作中深入研究。

参考文献：

[1] STROGATZ S H. Exploring complex networks[J]. Nature，2001，410（6825）：268-276.

[2] NéDA Z，RAVASZ E，VICSEK T，et al. Physics of the rhythmic applause[J]. Physical Review E，Statistical Physics，Plasmas，Fluids，and Related Interdisciplinary Topics，2000，61（6 Pt B）：6987-6992.

[3] NéDA Z，RAVASZ E，BRECHET Y，et al. The sound of many hands clapping. Tumultous applause can transform itself into waves of synchronized clapping[J]. Nature，2000，403（6772）：849-850.

[4] JOHN Buck E. Synchronous fireflies[J]. Scientific American，1976，234（5）：74-85.

[5] WIESENFELD K，COLET P，STROGATZ S H. Synchronization transitions in a disordered Josephson series array[J]. Physical Review Letters，1996，76（3）：404.

[6] STROGATZ S. Review of Sync：The emerging science of spontaneous order[J]. Notices of the American Mathematical Society，2003（3）.

[7] WIENER N. Nonlinear problems in random theory[J]. Students Quarterly Journal，1959，30（118）：77.

[8] WINFREE A T. Biological rhythms and the behavior of populations of coupled oscillators[J]. Journal of Theoretical Biology，1967，16（1）：15-42.

[9] KURAMOTO Y. Chemical Oscillations，Waves，and Turbulence[M]. Berlin，Heidelberg：Springer Berlin Heidelberg，1984：355-388.

[10] STROGATZ S H. From Kuramoto to Crawford：exploring the onset of synchronization in populations of coupled oscillators[J]. Physica D：Nonlinear Phenomena，2000，143（1/2/3/4）：1-20.

[11] KURAMOTO F，NISHIKAWA I. Onset of collective rhythms in large populations of coupled oscillators[C]//Kuramoto F，Nishikawa I. Springer Series in Synergetics. Berlin，Heidelberg：Springer Berlin Heidelberg，1989：300-306.

[12] TANAKA H A，LICHTENBERG A J，OISHI S. First order phase transition resulting from finite inertia in coupled oscillator systems[J]. Physical Review Letters，1997，78（11）：2104-2107.

[13] HONG H，CHOI M Y，YI J，et al. Inertia effects on periodic synchronization in a system of coupled oscillators[J]. Physical Review E，1999，59（1）：353-363.

[14] 黄霞，徐灿，孙玉庭，等. 耦合振子系统的多稳态同步分析[J]. 物理学报，2015，64（17）：170504-1-170504-11.

[15] VLASOV V，MACAU E E N，PIKOVSKY A. Synchronization of oscillators in a Kuramoto-type model with generic coupling[J]. Chaos：an Interdisciplinary Journal of Nonlinear Science，2014，24（2）：023120.

[16] KURAMOTO Y. International Symposium on Mathematical Problems in Theoretical Physics，Lecture Notes in Physics[C]//3rd IEEE International Symposium on Biomedical Imaging：Nano to Macro，2006 .

[17] CROSS M C，ZUMDIECK A，LIFSHITZ R，et al. Synchronization by nonlinear frequency pulling. [J]. Physical Review Letters，2004，93（22）：224101.

[18] MONTBRIó E，PAZó D. Shear diversity prevents collective synchronization[J]. Physical Review Letters，2011，106（25）：254101.

[19] CROSS M C，HOHENBERG P C. Pattern formation outside of equilibrium[J]. Review of Modern Physics，1993，65（3）：851-1112.

[20] STROGATZ S H，MIROLLO R E. Stability of incoherence in a population of coupled oscillators[J]. Journal of Statistical Physics，1991，63（3/4）：613-635.

[21] OTT E，ANTONSEN T M. Long time evolution of phase oscillator systems[J]. Chaos：an Interdisciplinary Journal of Nonlinear Science，2009，19（2）：023117.

[22] OTT E，ANTONSEN T M. Low dimensional behavior of large systems of globally coupled oscillators[J]. Chaos：an Interdisciplinary Journal of Nonlinear Science，2008，18（3）：037113.

[23] 陸启韶，彭临平，杨卓琴. 常微分方程与动力系统[M]. 北京：北京航空航天大学出版社，2010.

[责任编辑 杨 屹]