第３０卷第１２期 文章编号：１００６—９３４８（２０１３）１２—０２１９—０５ 计算机仿真 ２０１３年１２月 基于相互作用的耦合振子群体动力系统仿真 徐桂兰，唐宋，李清都 （重庆邮电大学非线性电路与系统研究所，重庆４０００６５） 摘要：大规模振子经耦合所形成的网络广泛存在于自然界的众多领域。由于存在高维非线性系统，研究网络的群体行为，是 对其进行设计控制和应用的前提。但是，对耦合振子网络群体构造的特定高维系统的研究往往涉及的大规模数值计算，算 法如果仅基于ＣＰＵ的计算平台往往会存在效率低的问题。为解决上述问题，提出了一种基于ＣＰＵ＋ＧＰＵ平台的异构并行 算法。算法首先设计了高效ＧＰＵ内核程序，然后借助ＧＰＵ加速，相对基于仅ＣＰＵ的平台，效率提高了近３０倍，实现了耦合 振子网络的高效仿真。 关键词：耦合振子网络；异构平台；并行计算 中图分类号：ＴＰ３０１．６ 文献标识码：Ｂ Ｓｉｍｕｌａｔｉｏｎ ｆｏｒ Ｄｙｎａｍｉｃ Ｓｙｓｔｅｍ Ｂａｓｅｄ ｏｎ Ｉｎｔｅｒａｃｔｉｎｇ Ｃｏｕｐｌｅｄ——Ｏｓｃｉｌｌａｔｏｒ Ｅｎｓｅｍｂｌｅｓ ＸＵ Ｇｕｉ—ｌａｎ，ＴＡＮＧ Ｓｏｎｇ，ＬＩ Ｑｉｎｇ—ｄｕ （Ｉｎｓｔｉｔｕｔｅ ｆｏｒ Ｎｏｎｌｉｎｅａｒ Ｃｉｒｃｕｉｔ ａｎｄ Ｓｙｓｔｅｍｓ，Ｃｈｏｎｇｑｉｎｇ Ｕｎｉｖｅｒｓｉｔｙ ｏｆ Ｐｏｓｔｓ ａｎｄ Ｔｅｌｅｃｏｍｍｕｎｉｃａｔｉｏｎｓ，Ｃｈｏｎｇｑｉｎｇ ４０００６５，Ｃｈｉｎａ） ＡＢＳＴＲＡＣＴ：Ｉｔ ｉｓ ａ ｗｉｄｅ ｒａｎｇｅ ｏｆ ｅｘｉｓｔｅｎｃｅ ｏｆ ｎｅｔｗｏｒｋｓ ｆｏｒｍｅｄ ｂｙ ｃｏｕｐｌｅｄ ｌａｒｇｅ—ｓｃａｌｅ ｏｓｃｉｌｌａｔｏｒｓ ｉｎ ｔｈｅ ｎａｔｕｒｅ． Ｔｈｅ ｓｔｕｄｙ ｏｆ ｔｈｅ ｇｒｏｕｐ ｂｅｈａｖｉｏｒ ｏｆ ｔｈｅ ｎｅｔｗｏｒｋｓ ｉｓ ａ ｐｒｅｍｉｓｅ ｔｈａｔ ｃａｒｒｉｅｓ ｏｎ ｄｅｓｉｇｎｉｎｇ，ｃｏｎｔｒｏｌｌｉｎｇ ａｎｄ ｐｕｔｔｉｎｇ ｉｎｔｏ ａｐ— ｐｌｉｃａｔｉｏｎ．Ｈｏｗｅｖｅｒ，ｔｈｅ ｐｒｏｂｌｅｍ ｏｆ ｍａｓｓｉｖｅ ｎｕｍｅｒｉｃａｌ ｃａｌｃｕｌａｔｉｏｎ ａｎｄ ｔｈｅ ｌｏｗ ｅｆｉｃｉｅｎｃｙ ｏｆ ｔｆｈｅ ｃｏｍｐｕｔｉｎｇ ｐｌａｔｆｏｒｍ ｂａｓｅｄ ｏｎ ＣＰＵ ａｙｅ ｏｆｔｅｎ ｉｎｖｏｌｖｅｄ ｉｎ ｔｈｅ ｒｅｓｅａｒｃｈ ｏｆ ｔｈｅ ｈｉｇｈ—ｄｉｍｅｎｓｉｏｎ ｓｙｓｔｅｍ ｂｕｉｌｔ ｂｙ ｔｈｅ ｃｏｕｐｌｅｄ ｏｓｃｉｌｌａｔｏｒ ｎｅｔｗｏｒｋ ｇｒｏｕｐｓ。Ｉｎ ｏｒｄｅｒ ｔｏ ｓｏｌｖｅ ｔｈｅ ｐｒｏｂｌｅｍｓ ｍｅｎｔｉｏｎｅｄ ａｂｏｖｅ，ｔｈｅ ｐａｐｅｒ ｆｉｒｓｔ ｐｒｅｓｅｎｔｅｄ ａ ｈｅｔｅｒｏｇｅｎｅｏｕｓ ｐａｒａｌｌｅｌ ａｌｇｏｒｉｓｍ ｂａｓｅｄ ｏｎ ＣＰＵ＆ＧＰＵ ｐｌａｔｆｏｒｍ ａｎｄ ｄｅｖｅｌｏｐｅｄ ａ ｈｉｇｈ ｅｆｉｃｉｅｎｔｆ ＧＰＵ ｋｅｒｎｅｌ ｐｒｏｇｒａｍ．Ｔｈｅｎ ｂｙ ａｄｄｉｎｇ ＧＰＵ ｉｎｔｏ ｔｈｅ ｐｌａｔ— ｆｏｒｍ，ｔｈｅ ｅｆｆｉｃｉｅｎｃｙ ｏｆ ｔｈｅ ｎｅｗ ｐｌａｔｆｏｒｍ ｗａｓ ｉｎｃｒｅａｓｅｄ ｎｅａｒｌｙ ３０ ｔｉｍｅｓ ｃｏｍｐａｒｅ ｔｏ ｔｈｅ ＣＰＵ ｐｌａｔｆｏｒｍ．Ｆｉｎａｌｌｙ，ｔｈｅ ｈｉｇｈ ｅｆｉｃｉｅｎｔｆ ｓｉｍｕｌａｔｉｏｎ ｏｆ ｔｈｅ ｃｏｕｐｌｅｄ ｏｓｃｉｌｌａｔｏｒ ｎｅｔｗｏｒｋ ｗａｓ ｉｍｐｌｅｍｅｎｔｅｄ． ＫＥＹＷＯＲＤＳ：Ｃｏｕｐｌｅｄ ｏｓｃｉｌｌａｔｏｒ ｎｅｔｗｏｒｋ；Ｈｅｔｅｒｏｇｅｎｅｏｕｓ ｐｌａｔｆｏｒｍ；Ｐａｒａｌｌｅｌ ｃｏｍｐｕｔｉｎｇ 群体行为，发挥着越来越重要的作用。例如：大群萤火虫的 １ 引言 随着非线性科学的不断发展，高维非线性系统的复杂行 同步闪烁…，大量观众的同频鼓掌 ，伦敦千禧桥的摇 摆 ｊ，细胞的无氧乳酸振荡 ，神经元协同震荡引起的癫 为逐渐成为研究的热点。该系统涉及很多领域，如流体动力 学、激光物理学、非线性光学、电子学和神经动力学…等等。 此类系统通常可以看做是大量相对简单元素的集体行为。 振子常被选为此类简单元素，因其之间具备不同类型的相互 痫 ］，电化学中大量振子的凝聚现象 等。研究这些群体行 为是对其设计控制和进行应用的前提。例如，制造高能效半 导体激光器 ］、基于Ｊｏｓｅｐｈｓｏｎ结阵列的大功率微波器件 和治疗癫痫 等。目前，国外关于振子网络群体行为的研究 较多，早期有日本的Ｋｕｒａｍｏｔｏ研究全局耦合振子 ，然后有 作用 ｊ，在许多情况下，系统的动力学行为可以在适当的低 维集体模式下得到有效体现。大规模的振子经简单耦合形 成的网络广泛存在于自然界的众多领域。网络所表现出的 基金项目：国家自然科学基金资助项目（６１１０４１５０，１０９７２０８２） 收稿日期：２０１３—０３—１５ Ｋｉｓｓ、Ｏｌｍｉ、ＭｃＣｕｌｌｅｎ、Ｓｔｒｏｇａｔｚ等也对其它振子网络群体进行 了深入研究 ’ ｔ”］，国内主要有王文旭、杨会杰等人研究Ｅ卜 ｄｏｓ—Ｒｅｎｙｉ网络 。 由于在实际情况中，对于耦合振子网络的研究一般存在 三点困难。第一，构成的特定高维系统的耦合振子网络群体 一２１９— 通常由成千上万个振子构成，其理论研究只能基于简化系 统，比如著名Ｋｕｒａｍｏｔｏ模型 ，但是，对于大多实际系统，这 种近似系统可能根本就不存在。因此，迫切需要开辟一些通 用方法来实施研究。第二，耦合振子网络一般是高维甚至是 超高维的系统，目前，并没有成熟的理论和研究方法。因此， 如何采用数值计算对群体混沌进行深入研究，成了一个极具 挑战的重要问题。第三，耦合振子网络个体数目的庞大导致 对其进行研究一般都会涉及大规模的数值计算，因此仅用 －（Ｑ ） 其中，ｃ和ｓ分别描述振子网络群体内部和相互之问的耦合， （ｔ）和 （ｔ）描述两个群体内部交替地进行着全局耦合。由 （５）式降阶后的控制方程为： ，ＣＰＵ做计算往往会出现速度慢，效率低和分辨率低等问题。 当今计算机构架面向并行计算和异构计算革新的趋势越来 警）＝ 越明显 。ＣＰＵ从单核发展到多核的过程中，多线程编程 ）＝（Ｑ ）鲁 ２ （ｆ）警＋ 也变得越来越重要。中央处理器ＣＰＵ和图形处理器ＧＰＵ都 有计算能力，前者拥有高效的串行指令，而后者有突出的并 行计算能力 １７ ３。异构计算则大多采用ＣＰＵ＋ＧＰＵ异构并 行。目前，随着ＧＰＵ的计算能力不断提高和快速发展，特别 是ＧＰＵ异构计算协议的制定 ，使得个人计算变得简单而 经济，且实现了计算机整体计算能力利用的最大化。 本文针对这些问题，提出了基于ＣＰＵ＋ＧＰＵ平台上的异 构算法。ＯｐｅｎＣＬ是一种可横跨ＣＰＵ、ＧＰＵ和其它处理器，让 软件开发人员的并行编程便携，且能高效地获取这些异构处 理平台能力的编程语言，因此本文采用ＭＡＴＬＡＢ、Ｃ＋＋和 ＯｐｅｎＣＬ相互协调的异构计算。为了更好地进行理论描述， 选用了 个简单而且方便的理想广义模型 。这个模型由 两个“ 的自维持振子网络群组成，当一个振子网络群体内 部进行全局耦合时，另一个振子群体内部的全局耦合关闭。＋ ＋ 振子群体之间附加的耦合是采用二次非线性平均场实现。甜 ｙ 然后对由此振子网络群体构造的特定高维系统在集体变量 ＝ ＝ 层面上进行仿真，最终证明其具有混沌行为。 ２耦合振子群体构造的高维系统的模型 ２．１＋ 一 基本广义模型 ＋ 一 考虑两个之间有相互作用，由自维持振子构造成的振子 网络群体。群体的大小都为Ｋ＝６４ ２ ｎ，群体中的成员分 ．吼 ．一 ｎ ｎ 别用 、ｆＹ 表示，其中ｋ的取值范围为１至Ｋ。在网络群体 内，振子之间通过平均场 和ｙ耦合， 和ｙ的表达式为： Ｘ＝ １ Ｙ＇＝－－…Ｚ ｙ （１） （１）式的导数ｄＸ／ｄｔ、ｄＹ／ｄｔ的表达式为： ｄＸ１ ｄＹ＝ ， ＝ （２） 网络群体问之问的非线性相互作用通过二阶平均场来 描述，表达式为： Ｘ２＝ １ ２ ＝ ｙ２ （３） （３）式的导数ｄＸ２／ｄｔ、ｄ ／ｄｆ的表达式为： ｄＸ２＝ ２ 鲁， ＝ ２ｙｋ ｄｙｋ ｃ４， 整个模型用以下二阶微分方程作为控制方程： ．．．．——２２０．．．．—— ｓｉｎ ０ｔ Ｊｙ（Ｙｋ， 、） ， ） ２ ＫＬ（ｔ） ＋ ｓｉｎ ∞ ｔ （６） ２．２模型的数值设置 为了对模型进行仿真以研究其性质，需要对广义模型 （５）中参数值进行具体设置。其中模型中振子的自然频率 ｔＯ 的分布如下： ＝ｔＯｏ—Ａ ＋ ａ盯 ｒｅｃｏｓ［Ｌ 凡 一 ］Ｊ （７） 对于（７）式，选择振子基频 。＝２１ｒ，进一步，令△ ＝ａｖ／８ 对其进行离散化。同时设置模型（５）中的 （ｔ）＝ＣＯＳ （竹∥ ） （ｔ）＝ｓｉｎ （￣ｒｔ／Ｔ），同步周期Ｔ＝１００，描绘振子群体内 部耦合的参数 ＝１以及描绘振子群体间耦合的参数 ＝ ０．１，单个振子幅度Ｑ＝３。 根据数值设置后的具体模型，设计了以下基于ＣＰＵ＋ ＧＰＵ异构计算的算法以对模型所构造的高维系统进行仿真， 算法的具体细节如下。 ３异构计算的系统仿真算法 ３．１算法的基本思想 ＭＡＴＬＡＢ具备强大的作图功能，本算法基本思想是通过 建立ＭＥＸ文件，使ＭＡＴＬＡＢ能够调用Ｃ＋＋子程序，从而分 别发挥ＭＡＴＬＡＢ和ＧＰＵ的各自优势。首先在ＭＡＴＬＡＢ中进 行系统初始化，令Ｋ＝６４ ２ ｎ，设置绘图的最大时问 …： ２００、采样间隔（同时也是时间间隔）ｈ＝０．０１。计算出采样点 数后，将其赋值给Ｍ。将系统参数，绘图参数以及数１／６，１／ ３，１／２等赋给一个３０ ｘ １的数组Ｐ。系统初始化即构造振子 、Ｙ 及其导数ｄｘ ／ｄｔ、ｄｙ＾／ｄｔ的初始值且赋值给Ｋ×４的矩 阵 。然后调用ＭＥＸ函数和处理其接口，将由系统初始值， 频率值所组成的数组以及Ｐ作为ＭＥＸ函数的输入值。把第 一个输入变量，即系统初始值和频率值所组成的数组的首地 址赋给ｐｚ；把第二输入变量数组Ｐ的首地址赋给数组ＰＰ；同 时创建一个Ｍ×５的矩阵，用以存放输出数据，同时将矩阵首 地址赋给ＰＺ。最后进入ＯｐｅｎＣＬ宿主程序部分。宿主程序 共分为三部分：第一部分是对ＯｐｅｎＣＬ进行初始化，包括计算 平台（Ｐｌａｔｆｏｒｍ）的寻找，上下文（Ｃｏｎｔｅｘｔ）的创建，对计算设备 ３）将ｔ移动半个步长，程序转人利用计算出的第二阶的 点求（２）和（４）式。然后进人第二个－厂０ｒ循环求Ｒ—Ｋ法第二 阶相应值。 。和 ：用与第一阶相同的方法计算。 ［ｉ］和 ５ ［ｉ］在第一阶值的基础上分别叠加上１／３乘以 和 ， 同时将（６）式第三阶要使用到的点求出来。第三个 ｒ循环 （Ｄｅｖｉｃｅｓ）的探测和创造ＯｐｅｎＣＬ命令队列（Ｃｏｍｍａｎｄ Ｑｕｅｕｅ）。第二部分是建立ＧＰＵ内核对象，如果ＧＰＵ的内核 程序是二进制文件，则直接加载，然后创建一个程序对象。 若不是二进制文件而是源程序，则创建程序对象，用该程序 对象创建内核对象，同时输出二进制内核文件，以便下次调 按照与第二个循环相似的方法计算 、ｒｋ 、ＳＶ［ｉ］、 ［ｉ］。 因为在计算（６）式第四阶将要使用到的点时只需用第三阶斜 率乘以ｈ，所以算法中仅用当次最外层 ｒ循环的初始值加上 ｌ用。第三部分是运行内核程序，包括获取计算设备的具体信 息、得到给定名称的内核对象句柄以及数据传输与并行计 和 ２。 ４）时间ｔ移动半个步长进入第四个 ｒ循环。计算方法 算。然后执行内核，完成后将数据读回，并传回ＭＡＴＬＡＢ，再 与上述类似。根据Ｒ—Ｋ法，叠加到ＳＶ［ｉ］、ＳＶ ［ｉ］的值用１／６ 在ＭＡＴＬＡＢ中进行曲线绘制。 整个算法的核心部分，也是体现ＧＰＵ并行算法的地方 是内核程序，下面详细介绍内核程序。 ３．２算法的结构与内容 ３．２．１程序初始化 ＯｐｅｎＣＬ宿主程序中，为了让常数便于存取，在ＧＰＵ的 常数内存区创建只读缓存区ｚｂ和ｐｂ，在ＧＰＵ的全局内存区 创建可读可写缓存区ｚＢ，并将指针ｐｚ，ｐｐ，ＰＺ所指向的内容 分别复制给指针ｚｂ，ｐｂ和ｚＢ．然后将它们设置为内核对象的 第一、第二和第三参数。同时将线程总数及每一个线程组大 小均设置为１２８。然后进入内核程序部分。首先对（６）进行 宏定义，在ＧＰＵ的本地内存区创建缓存区ｂｆ和ｓｂｆ。对该程 序将要用到的变量进行初始定义后，便将数组ｐｂ中的参数 值赋给相应变量以便以后调用，此时，数据ｐ６为二维空间中 的数组。读取参数后，利用线程的，ｚ）且一个线程处理ｎ个 数据，根据 。等于ｚ６加上线程 乘以ｎ所体现的并行计 算，后面的很多数据处理都只需要考虑一个线程的ｎ个振子 即可． ３．２．２龙格库塔算法部分 整个内核程序的最终目的是用一个外层 ｒ循环加四个 内层 ｒ循环实现系统（６）的龙格库塔算法（以下简称Ｒ—Ｋ 法）。详细流程如下： １）外层循环共ｚ层。首先根据传递的第Ｚ×ｈ时刻的系 统值计算（１）、（２）和（４）式（计算出（４）式可使（６）式中的 （‰，ｄｘ ／ｄｔ）√ （Ｙ ，ｄｙ ／ｄｔ）相互独立，便于后面计算 —Ｋ 法中的四个斜率）。其中ｚ为循环变量且初值为Ｚ＝０。将 计算出的Ｘ、ｄＸ／ｄｔ、Ｙ、ｄＹ／ｄｔ值分别存在Ｐｚ的第ｚ＋１行。 ２）进入内层，０ｒ循环。每个内层循环为ｎ层。第一个 －厂０ｒ循环用来Ｒ—Ｋ法第一阶相应的值：因计算第二阶斜率时 要使用第一阶斜率乘以ｈ／２，且计算系统的下一个点也要乘 以矗，为了简便起见，算法提前乘以 。这些相应的值用 和 ２表示， ，的两维分别等于 （ ，ｄｘ ／ｄｔ）￣ｆｙ（Ｙ ， ／ ｄｔ）乘以ｈ，同理 的两维等于厶（ ， ／ ）和厶（Ｙ ， ／ｄｔ）乘以ｈ。然后根据Ｒ—Ｋ法计算系统下一个点的公 式，将 ， ：乘以１／６后，用定义的ＳＶ［ｉ］、ＳＶ ｉ］存起来，同 时将求Ｒ—Ｋ法第二阶相应值将要使用到的点计算出来。 乘以 ．和 ，再叠加上当次外层 ｒ循环的初始值以得到 系统的间隔为ｈ的下一个点。将它们分别赋给 ［ｉ］和 ［ｉ］ 来作为求下一次循环的初始值。 ５）最后利用计算所得的下一次循环的初始值平均场 、 ｌ，及其导数ｄＸ／ｄｔ、ｄＹ／ｄｔ的下一个点。将其与相应的时间 （时间间隔为ｈ）分别存在Ｐｚ的第ｚ＋２行。然后Ｚ自加１， 直到完成Ｍ一１次循环。将得到的矩阵返值给宿主程序。 算法的广义详细控制流程见图１。 图１ Ｒ—Ｋ法的广义控制流程图 ３．２．３算法的仿真结果与性能分析 图２至５为当单个网络群体大小 取１２８，群体个数为 ２时平均场 、ｙ及其导数ｄＸ／ｄｔ、ｄＹ／ｄｔ最终的仿真结果。从 仿真图１和３可以看出，当调制 （ｔ）和 （ｔ）的值较小时， 平均场 、ｙ几乎为０，当 （￡）和 （ｆ）值较大时， 、ｌ，则变 为最大。此仿真结果与文献［１９］的仿真结果是完全一致的。 从图２和图４可以发现，ｄＸ／ｄｔ、ｄＹ／ｄｔ都稳定到了周期轨，即 最终都进入了稳态，说明网络群体的成员最后同步。 改变Ｋ＝６４ ２ ｎ中的ｎ，使之为１０，图６至９为当群 体个数为２０时 、ｙ及它们的导数的仿真结果。比较图２—５ 和图５—９，可见虽然修改了个数和提高了速度，最终图的结 果是一致的。 一２２ｌ一 图９平均场Ｙ的导数的仿真结果图 表１算法测试结果对比 ４结论 本文提出了基于ＣＰＵ＋ＧＰＵ异构平台，对耦合振子网络 群体构造的特定高维系统所涉及的大规模数值计算，相对于 一般的非异构算法，提高了接近３０倍以上的计算速度．并行 是把一个复杂问题，分解成多个能同时处理的子问题，非常 适合具有大规模点的快速计算．本文虽然只对大小为１２８、 １２８０和１５３６的振子网络群进行仿真，对于具有更多群体个 数的耦合振子网络，此方法可通过直接修改Ｋ的值、增加群 体个数和每个线程处理个数得以沿用。而且，基于ＣＰＵ＋ ＧＰＵ异构平台是计算机未来的发展方向，因此本文所涉及到 的算法的应用前景还是十分广阔，它能够解决了一些大规模 计算问题，如基于ＣＰＵ＋ＧＰＵ异构平台应用于吸引盆，可使 求出的吸引盆更加细致，包含更多信息，使对吸引盆的数值 研究变得更加细致全面。 参考文献： ［１］Ａ Ｐｉｋｏｖｓｋｙ，ｅｔ ａ１．Ｓｙｎｃｈｒｏｎｉｚａｔｉｏｎ：ａ ｕｎｉｖｅｒｓａｌ ｃｏｎｃｅｐｔ ｉｎ ｎｏｎｌｉｎ— ｅａｒ ｓｃｉｅｎｃｅ［Ｍ］．Ｃａｍｂｒｉｄｇｅ ｕｎｉｖｅｒｓｉｔｙ ｐｒｅｓｓ，２００３． ［２］ Ｚ Ｎｒｄａ，ｅｔ ａ１．Ｔｈｅ ｓｏｕｎｄ ｏｆ ｍａｎｙ ｈａｎｄｓ ｃｌａｐｐｉｎｇ—Ｔｕｍｕｌｔｕｏｕｓ 印ｐｌａｕｓｅ ｃａｎ ｔｒａｎｓｆｏｒｍ ｉｔｓｅｌｆ ｉｎｔｏ ｗ￣ｖｅｓ ｏｆ ｓｙｎｃｈｒｏｎｉｚｅｄ ｃｌａｐｐｉｎｇ ［Ｊ］．Ｎａｔｕｒｅ，２０００，４０３（６７７２）：８４９—８５０． ［３］ Ｓ Ｈ Ｓｔｒｏｇａｔｚ，ｅｔ ａ１．Ｔｈｅｏｒｅｔｉｃａｌ ｍｅｃｈａｎｉｃｓ：Ｃｒｏｗｄ ｓｙｎｃｈｒｏｎｙ ｏＮ ｔｈｅ Ｍｉｌｌｅｎｎｉｕｍ Ｂｒｉｄｇｅ［Ｊ］．Ｎａｔｕｒｅ，２００５，４３８（７０６４）：４３—４４． ［４］Ｄ Ｇｏｎｚｅ，Ｎ Ｍａｒｋａｄｉｅｕ，Ａ Ｇｏｌｄｂｅｔｅｒ．Ｓｅｌｅｃｔｉｏｎ ｏｆ ｉｎ—ｐｈａｓｅ ｏｒ ｏｕｔ—－ｏｆ——ｐｈａｓｅ ｓｙｎｃｈｒｏｎｉｚａｔｉｏｎ ｉｎ ａ ｍｏｄｅｌ ｂａｓｅｄ ｏｎ ｇｌｏｂａｌ ＣＯＵ－ ｐｌｉｎｇ ｏｆ ｃｅｌｌｓ ｕｎｄｅｒｇｏｉｎｇ ｍｅｔａｂｏｌｉｃ ｏｓｃｉｌｌａｔｉｏｎｓ［Ｊ］．Ｃｈａｏｓ：Ａｎ Ｉｎ— ｔｅｒｄｉｓｅｉｐｌｉｎａｒｙ Ｊｏｕｒｎａｌ ｏｆ Ｎｏｎｌｉｎｅａｒ Ｓｃｉｅｎｃｅ，２００８，１８（３）：０３７１２７ —０３７１２７—１２． ［５］ Ｉ Ｏｓｏｒｉｏ．ｅｔ ａ１．Ｅｐｉｌｅｐｓｙ：Ｔｈｅ Ｉｎｔｅｒｓｅｃｔｉｏｎ ｏｆ Ｎｅｕｒｏｓｃｉｅｎｃｅｓ，Ｂｉｏｌ— ｏｇｙ，Ｍａｔｈｅｍａｔｉｃｓ，Ｅｎｇｉｎｅｅｒｉｎｇ ａｎｄ Ｐｈｙｓｉｃｓ［Ｍ］． ＣＲＣ Ｐｒｅｓｓ，２０１１． ［６］ Ｉ Ｚ Ｋｉｓｓ．Ｙ Ｚｈａｉ．Ｊ Ｌ Ｈｕｄｓｏｎ．Ｅｍｅｒｇｉｎｇ ｃｏｈｅｒｅｎｃｅ ｉｎ ａ ｐｏｐｕｌａ－ ｔｉｏｎ ｏｆ ｃｈｅｍｉｃａｌ ｏｓｃｉｌｌａｔｏｒｓ［Ｊ］．Ｓｃｉｅｎｃｅ，２００２，２９６（５５７３）：１６７６ １６７８． ［７］ Ｍ Ｃ Ｃｒｏｓｓ．Ｐ Ｃ Ｈｏｈｅｎｂｅｒｇ．Ｐａｔｔｅｒｎ ｆｏｒｍａｔｉｏｎ ｏｕｔｓｉｄｅ ｏｆ ｅｑｕｉｌｉｂｒｉ— ａｍ［Ｊ］．Ｒｅｖｉｅｗｓ ｏｆ Ｍｏｄｅｍ Ｐｈｙｓｉｃｓ，１９９３，６５（３）：８５１． ［８］ Ｙ Ｋｕｒａｍｏｔｏ．Ｓｅｌｆ—ｅｎｔｒａｉｎｍｅｎｔ ｏｆ ａ ｐｏｐｕｌａｔｉｏｎ ｏｆ ｃｏｕｐｌｅｄ ｎｏｎ— ｌｉｎｅａｒ ｏｓｃｉｌｌａｔｏｒｓ［ｃ］．Ｉｎｔｅｒｎａｔｉｏｎａｌ ｓｙｍｐｏｓｉｕｍ ｏｎ ｍａｔｈｅｍａｔｉｃａｌ ｐｒｏｂｌｅｍｓ ｉｎ ｔｈｅｏｒｅｔｉｃａｌ ｐｈｙｓｉｃｓ：Ｓｐｒｉｎｇｅｒ，１９７５：４２０—４２２． ［９］ Ｍ Ｎｉｘｏｎ，ｅｔ ａ１．Ｓｙｎｃｈｒｏｎｉｚｅｄ ｃｌｕｓｔｅｒ ｆｏｒｍａｔｉｏｎ ｉｎ ｃｏｕｐｌｅｄ ｌａｓｅｒ ｎｅｔｗｏｒｋｓ［Ｊ］．Ｐｈｙｓｉｃａｌ Ｒｅｖｉｅｗ Ｌｅｔｔｅｒｓ，２０１１，１０６（２２）：２２３９０１． ［１Ｏ］ Ｍ Ｂｅｎｎｅｔｔ，Ｋ Ｗｉｅｓｅｎｆｅｌｄ．Ａｖｅｒａｇｅｄ ｅｑｕａｔｉｏｎｓ ｆｏｒ ｄｉｓｔｉｒｂｕｔｅｄ Ｊｏ－ ｓｅｐｈｓｏｎ ｊｕｎｃｔｉｏｎ ａｒｒａｙｓ［Ｊ］．Ｐｈｙｓｉｃａ Ｄ：Ｎｏｎｌｉｎｅａｒ Ｐｈｅｎｏｍｅｎａ， ２００４，１９２（３）：１９６—２１４． ［１１］ Ｓ Ｈ Ｓｔ Ｓｔｒｏｇａｔｚ．Ｆｒｏｍ Ｋｕｒａｍｏｔｏ ｔｏ Ｃｒａｗｆｏｒｄ：ｅｘｐｌｏｒｉｎｇ ｔｈｅ ｏｎｓｅｔ ｏｆ ｓｙｎｃｈｒｏｎｉｚａｔｉｏｎ ｉｎ ｐｏｐｕｌａｔｉｏｎｓ ｏｆ ｃｏｕｐｌｅｄ ｏｓｃｉｌｌａｔｏｒｓ［Ｊ］．Ｐｈｙｓｉ— ｃａ Ｄ：Ｎｏｎｌｉｎｅａｒ Ｐｈｅｎｏｍｅｎａ，２０００，１４３（１）：ｌ一２０． ［１２］ Ｎ ＭｃＣｕｌｌｅｎ，Ｐ Ｍｏｒｅｓｃｏ．Ｒｏｕｔｅ ｔｏ ｈｙｐｅｒｃｈａｏｓ ｉｎ ａ ｓｙｓｔｅｍ ｏｆ ＣＯＨ— ｐｌｅｄ ｏｓｃｉｌｌａｔｏｒｓｗｉｔｈｍｕｈｉｓｔａｂｉｌｉｔｙ［Ｊ］．ＰｈｙｓｉｃａｌＲｅｖｉｅｗ Ｅ，２０１１， ８３（４）：０４６２１２． ［１３］ Ｓ Ｏｌｍｉ，Ａ Ｐｏｌｉｔｉ，Ａ Ｔｏｒｃｉｎｉ．Ｃｏｌｌｅｃｔｉｖｅ ｃｈａｏｓ ｉｎ ｐｕｌｓｅ—ｃｏｕｐｌｅｄ ｎｅｕｒａｌ ｎｅｔｗｏｒｋｓ［Ｊ］．ＥＰＬ（Ｅｕｒｏｐｈｙｓｉｃｓ Ｌｅｔｔｅｒｓ），２０１１，９２ （６）：６０００７． ［１４］ Ｈ Ｙａｎｇ，Ｆ Ｚｈａｏ，Ｂ Ｗａｎｇ．Ｃｏｌｌｅｃｔｉｖｅ ｃｈａｏｓ ｉｎｄｕｃｅｄ ｂｙ ｓｔｒｕｅ－ ｔｕｒｅｓ ｏｆ ｃｏｍｐｌｅｘ ｎｅｔｗｏｒｋｓ［Ｊ］．Ｐｈｙｓｉｃａ Ａ：Ｓｔａｔｉｓｔｉｃａｌ Ｍｅｃｈａｎｉｃｓ ａｎｄ ｉｔｓ Ａｐｐｌｉｃａｔｉｏｎｓ，２００６，３６４：５４４—５５６． ［１５ ］Ｊ Ａ Ａｃｅｂｒ６ｎ，ｅｔ ａ１．Ｔｈｅ Ｋｕｒａｍｏｔｏ ｍｏｄｅｌ：Ａ ｓｉｍｐｌｅ ｐａｒａｄｉｇｍ ｏｆｒ ｓｙｎｃｈｒｏｎｉｚａｔｉｏｎ ｐｈｅｎｏｍｅｎａ［Ｊ］．Ｒｅｖｉｅｗｓ ｏｆ ｍｏｄｅｍ ｐｈｙｓｉｃｓ， ２００５，７７（１）：１３７． ［１６］ Ｗ Ｈｗｕ，Ｄ Ｋｉｒｋ．Ｐｒｏｇｒａｍｍｉｎｇ ｍａｓｓｉｖｅｌｙ ｐａｒｌａｌｅｌ ｐｒｏｃｅｓｓｏｒｓ［Ｊ］． Ｓｐｅｃｉａｌ Ｅｄｉｔｉｏｎ，２００９，９２． ［１７］ 李清都，周红伟，杨晓松．基于异构计算的简单行走模型的吸 引区域研究［Ｊ］．物理学报，２０１２，６１（４）：４０５０３—０４０５０３． ［１８］ Ｇｒｏｕｐ Ｋｏｗ．Ｔｈｅ ＯｐｅｎＣＬ Ｓｐｅｃｉｉｆｃａｔｉｏｎ Ｖｅｒｓｉｏｎ １．１，２０１０［Ｊ］． Ｄｏｃｕｍｅｎｔ Ｒｅｖｉｓｉｏｎ，２０１０，４４． ［１９］ Ｓ Ｐ Ｋｕｚｎｅｔｓｏｖ，Ａ Ｐｉｋｏｖｓｋｙ，Ｍ Ｒｏｓｅｎｂｌｕｍ．Ｃｏｌｌｅｃｔｉｖｅ ｐｈａｓｅ ｃｈａ— ｏｓ ｉｎ ｔｈｅ ｄｙｎａｍｉｃｓ ｏｆ ｉｎｔｅｒａｃｔｉｎｇ ｏｓｃｉｌｌａｔｏｒ ｅｎｓｅｍｂｌｅｓ［Ｊ］．ａｒＸｉｖ ｐｒｅｐｒｉｎｔ ａｒＸｉｖ：１０１２０６５２，２０１０． ［作者简介］ 徐桂兰（１９８５一），女（汉族），四ＪＩｌ省遂宁市人，硕 士研究生，主要研究领域为动力系统、数值计算。 唐宋（１９８２一），男（汉族），四川省达州人，硕士 研究生，主要研究领域为动力系统、数值计算。 李清都（１９８０一），男（汉族），重庆人，教授，博士， 主要研究方向：混沌动力系统、流形计算、拓扑马蹄理论。 －－－－——２２３．．．——